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Abstract 

The bottom quark pole mass Mb is determined using a sum rule which relates the masses and 
the electronic decay widths of the T mesons to large n moments of the vacuum polarization 
function calculated from nonrelativistic quantum chromodynamics. The complete set of ncxt-to- 
next-to-leading order (i.e. 0(a 2 s , a s v, v 2 ) where v is the bottom quark cm. velocity) corrections 



is calculated and leads to a considerable reduction of theoretical uncertainties compared to a pure 
next-to-leading order analysis. However, the theoretical uncertainties remain much larger than the 
experimental ones. For a two parameter fit for Mb, and the strong MS coupling a s , and using the 
scanning method to estimate theoretical uncertainties, the next-to-next-to-leading order analysis 
yields 4.74 GeV < M b < 4.87 GeV and 0.096 < a s {M z ) < 0.124 if experimental uncertainties 
are included at the 95% confidence level and if two-loop running for a s is employed. Mb and 
a s have a sizeable positive correlation. For the running MS bottom quark mass this leads to 
4.09 GeV < mf,(Mx(is)/2) < 4.32 GeV. If a s is taken as an input, the result for the bottom 
quark pole mass reads 4.78 GeV < M b < 4.98 GeV (4.08 GeV < m b (M T(ls) /2) < 4.28 GeV) for 
0.114 < a s (M z ) < 0.122. The discrepancies between the results of three previous analyses on the 
same subject by Voloshin, Jamin and Pich, and Kuhn et al. are clarified. A comprehensive review 
on the calculation of the heavy quark-antiquark pair production cross section through a vector 
current at next-to-next-to leading order in the nonrelativistic expansion is presented. 



1 Introduction 



Quantum chromo dynamics (QCD) is the established theory of the strong interactions. The determina- 
tion of its parameters, the strong coupling and the quark masses, and continuous tests of its consistency 
with experimental measurements belong to the most important tasks within particle physics. For the 
strong coupling an almost countless number of determinations exists. The most precise determinations 
now quote uncertainties in a s (M z ) of less than 5%.[| The remarkable feature of the a s determinations, 
however, is their consistency to each other (see e.g. for a review). The situation for the quark masses 
can certainly be described as much less coherent. For the bottom quark pole mass, which represents 
an important ingredient for the theoretical description of B mesons decays the determination of the 
corresponding Cabibbo-Kobayashi-Maskawa matrix elements, the situation is particularly confusing. 
In the past few years there have been three determinations by Voloshin (Mj, = 4.827 ± 0.007 GeV) 
and later by Jamin and Pich (M b = 4.60 ± 0.02 GeV)|] and Kiihn et al. (M b = 4.75 ± 0.04 GeV) @ 
which, although they have all been obtained from the same experimental data on the spectrum and 
the electronic decay widths of the T mesons, are contradictory to each other if the quoted uncertain- 
ties are taken seriously. Further, the three analyses |2|, |3|, |j were all based on the same sum rule 
which relates large n moments (i.e. large number of derivatives at zero momentum transfer) of the 
vacuum correlator of two bottom-antibottom vector currents to an integral over the total production 
cross section of hadrons containing a bottom and an antibottom quark in e + e~ annihilation. In the 
limit of large n the moments can be calculated in a nonrelativistic expansion || ^] and higher order 
(relativistic) corrections can be implemented in a systematic way. 

This paper contains a determination of the bottom quark pole mass, where theoretical un- 
certainties are treated in a conservative way. It is partly motivated by the belief that a carefully 
performed analysis of theoretical uncertainties is mandatory in order to see whether the uncertainties 
presented in Refs. |^, ||, Q are realistic. In this work the method of choice is to scan all theoretical 
parameters independently over reasonably large windows. We will show that this method to estimate 
theoretical uncertainties is more conservative than the methods used in Refs. ||, [3|, In particular it 
renders the results obtained by Voloshin and Kiihn et al. consistent to each other. With the scanning 
method the precise results of Voloshin and Kiihn et al. can only be obtained if some model-like as- 
sumptions are imposed which are beyond first principles QCD. The by far bigger part of the motivation 
for this work, however, comes from the fact that now the technical and conceptual tools have been 
developed 0, || |j| [l(J to include the next-to-next-to-leading order (NNLO) relativistic corrections to 
the large n moments into the analysis. A large fraction of this paper is devoted to a comprehensive 
presentation and review of the concepts and calculations necessary to determine those NNLO contri- 
butions. In particular, we use the concept of effective field theories formulated in the frame work of 



nonrelativistic quantum chromodynamics (NRQCD) [Tl], 12] to deal with the problems of ultraviolet 
divergences which arise if relativistic corrections to the expressions in the nonrelativistic limit are 
calculated. However, we regard NRQCD merely as a technical tool and do not spend too much time 
on formal considerations. Whenever possible we rely on physical rather than formal arguments and 
use results from older literature even if they have not been derived in the framework of NRQCD. It 
is the main intention of this work to calculate the NNLO corrections to the large n moments and 
to analyze their impact on the determination of M b . We show that the NNLO corrections lead to a 



Throughout this paper the strong coupling is defined in the MS scheme. 
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considerable reduction of theoretical uncertainties in the determination of Mj,. 

The program of this paper is as follows: In Section |2| we introduce our notation and explain the 
ideas and concepts on which our analysis and calculations are based on. NRQCD is introduced and a 
recipe for the calculation of the moments at NNLO is presented. Because the heavy quark-antiquark 
cross section in the threshold regime represents an important intermediate step in the calculation 
of the moments, Section || also contains a comprehensive review on the basic concepts involved in 
the calculation of the vector current induces cross section at NNLO. In Section || all calculations are 
carried out explicitly and all relevant formulae are displayed. Section ||] contains a discussion on some 
peculiarities of the large n moments. A detailed description of the treatment of the experimental data, 
the fitting procedure and the scanning method is given in Section |5[ In Section ^ the numerical results 
are presented and discussed. Two different determinations of M\, are carried out. First, M& and a s 
are fitted simultaneously and, second, is fitted while a s is taken as in input. In Section 0, finally, 
we comment on the three previous analyses in Refs. [0, [3|, ||] and Section |8| contains the conclusions. 
Attached to this paper are three appendices which contain material which we found too detailed to 
be presented in the main body of the paper. 

The reader who is mainly interested in the results for the bottom quark mass can safely skip 
most of Section [|, and Sections || and |] completely. 



2 The Basic Ideas and Notation 

The Sum Rule 

We start our consideration from the correlator of two electromagnetic currents of bottom quarks at 
momentum transfer q 

IV(g) = -iJdxe^iOlTj^jt^lO), (1) 

where 

jj(x) =6(x) 7m 6(x). (2) 

The symbol b denotes the bottom quark Dirac field. We define the n'th moment P n of the vacuum 
polarization function as 

(3) 

t2 



p. = *4$(£) v(«) 



q 2 =0 



nlq 2 \dq 2 , 

where Qb = —1/3 is the electric charge of the bottom quark. Due to causality the n-th moment P n 
can be written in terms of a dispersion integration 

ds 



where 

R(s) = 2 °1± (5) 

is the total photon mediated cross section of bottom quark-antiquark production in e + e~ annihilation 
normalized to the point cross section a v t = Ana 2 /3s. Assuming global duality, P n can be either 
calculated from experimental data for the total cross section in e + e~ annihilation^ or theoretically 
2 At the level of precision in this work the Z mediated cross section can be safely neglected. 
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using quantum chromodynamics (QCD). It is the basic idea of this sum rule to set the moments 
calculated from experimental data, P^ 1 , equal to those determined theoretically from QCD, P* , and 
to use this relation to determine the bottom quark mass (and the strong coupling) by fitting theoretical 
and experimental moments for various values of n. || |(| 

At this point it is mandatory to discuss the range of n for which the theoretical moments can 
be calculated sufficiently accurate (using perturbative QCD) to allow for a reliable extraction of M b 
and a s . Prom Eq. (Q) is obvious that each moment P n effectively corresponds to a smearing of the 
cross section R over some energy region AE located around the threshold point. Thus, only if the 
smearing range is sufficiently larger than A QCD ~ 0(200 — 300 MeV), a perturbative calculation of 
the moments is feasible Jl4|]. [In Ref. [14| is was argued that AE should be larger than AM b a s to 
avoid the complications involving a resummation of the Coulomb singularities oc (a s /v) m . Because 
this resummation is explicitly carried out at the NNLO level in this work, we have to take A QCD , the 
typical hadronic scale, as the size of the minimal smearing range.] We therefore conclude that n is 
not allowed to be too large if perturbative QCD shall be employed. We can derive an approximate 
upper bound for the allowed values of n by changing the integration variable in relation (Q) to the 



energy E 
expression (0 

Pn 



2M b . For n 3> 1 only energies E <C M b contribute, which allows us to expand 



for small E/M b (while regarding {E/M b )n of order one) 



n>l 



(4M ft 2 ) 



dE 
M~ b 



exp 



E 



■ n 



R[(2M b + E)'- 



,„/ E E 2 
1 + — , — on 
\M b Ml 



From Eq. @ we see that the size of the smearing range AE for large n is or order M b /n, 



AE 



Mb 
n 



(6) 



(7) 



Demanding that AE is larger than A QCD yields that the values of n for which a perturbative calculation 
of the moments can trusted should be sufficiently smaller than 15 — 20. To avoid systematic theoretical 
errors as much as possible we take 

n-max = 10 (8) 

as the maximal value for n employed in this work. On the other hand, it is also desirable to choose 
n as large as possible because the experimental cross section for electron positron annihilation into bb 
hadrons is much better known in the T resonance regime y/s ~ 9.5 — 10.5 GeV than above the BB 
threshold. By taking n large the lower lying resonance contributions in Eq. (|j) are enhanced relative 
to the continuum contributions leading effectively to a suppression of the experimental uncertainties 
in the continuum cross section |(| . For our analysis we choose 



n. 



(9) 



as the minimal value for n. It is the regime 4 < n < 10 which we will refer to as "large n" in this 
work. It is a very important fact that for 4 < n < 10 the bottom-antibottom quark dynamics in the 
theoretical moments P^ h is already nonrelativistic in nature. This can be seen by once again examining 
relation (|6|). Because for a given value of n only energies E < M b /n contribute, the corresponding 
bottom quark velocities v = \JE/M b (in the cm. frame) are in the range \v\ <^ 0.5, i.e. they are 
always considerably smaller than the speed of light. In particular, the velocity is already as large as 
the typical size of the strong coupling a s (M b v) ~ 0.3 governing the exchange of longitudinal polarized 
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gluons (in Coulomb gauge) among the bottom-antibottom quark pair. This leads to the breakdown 
of the conventional multi-loop perturbation expansion because the exchange of m longitudinal gluons 
generates singular terms oc (a s /v ) m , m = 0, 1, 2, . . ., (Coulomb singularities) in the cross section for 
small velocities. These singular terms would have to be resummed^to all orders in multi-loop pertur- 
bation theory in order to arrive at a viable description of the bottom-antibottom quark dynamics. In 
other words, the Coulomb interaction between the bottom and the antibottom quark has to be treated 
exactly ||. Because this is an impossible talk in the framework of covariant multi-loop perturbation 
theory it is mandatory to calculate the cross section and the theoretical moments in the nonrelativistic 
approximation by solving the Schrodinger equation supplemented by relativistic corrections. 



Perturbative NRQCD and the Cross Section 



In this paper we use NRQCD Ej , 12] to set up a consistent framework in which the corrections to 
the nonrelativistic limit (in form of the nonrelativistic Schrodinger equation) can be determined in 
a systematic manner at NNLO. This corresponds to corrections up to order a?, a s v and v 2 to the 
expressions in the nonrelativistic limit. We count orders of a s as orders of v because we treat the bb 
system as Coulombic. In the framework of multi-loop perturbation theory this would correspond to 
a resummation of all terms oc a™v k with m + k = 1,2,3 in the cross section for the small velocity 
expansion. NRQCD is an effective field theory of QCD designed to handle nonrelativistic heavy-quark- 
antiquark systems to in principle arbitrary precision. NRQCD is based on the separation of long- and 
short-distance effects by reformulating QCD in terms of an unrenormalizable Lagrangian containing 
all possible operators in accordance to the symmetries in the nonrelativistic limit. Treating all quarks 
of the first and second generation as massless and taking into account only those terms relevant for 



the NNLO calculation in this work the NRQCD Lagrangian reads [12] 



£nrqcd = ~TtG^G^ + £ 9ipQ 



q=u,d,s,c 



D D 



l 



iD t + ai — — + a 2 



2M t 8M t 3 



2M t 8M f 2 V ' 8M ( 2 V ' 



tp + 



+XX^ bilinear terms and higher dimensional operators . (10) 

The gluons and massless quarks are described by the conventional relativistic Lagrangian, where 
G^v is the gluon field strength tensor, q the Dirac spinor of a massless quark and the gauge 
covariant derivative. For convenience, all color indices in Eq. ( JT0| ) and throughout this work are 
suppressed. The nonrelativistic bottom and antibottom quarks are described by the Pauli spinors ip 
and X) respectively. Dt and D are the time and space components of the gauge covariant derivative 
D and E l = G° l and B i = ^^ k G^ k the electric and magnetic components of the gluon field strength 
tensor (in Coulomb gauge). The straightforward x*X bilinear terms are omitted and can be obtained 
using charge symmetry. The short-distance coefficients ai, . . . , 05 are normalized to one at the Born 

3 In this context "resummation" technically means that one carries out the resummation of singular terms in the 
(formal) kinematic regime a 3 \v\. The resulting series (uniquely) define analytic functions which can then be continued 
to the regime \v\ ^ a s . 
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level. The actual form of the higher order contributions to the short-distance coefficients a\, . . . ,05 
(and also to b\, b 2 in Eq. (|f~2l)) is irrelevant for this work, because we will later use the "direct matching" 
procedure Jj], at the level of the final result for the cross section. 

Let us first discuss the cross section R in the nonrelativistic regime. To formulate R in the 
nonrelativistic regime at NNLO in NRQCD we start from the fully covariant expression for the total 
cross section 



R(i 2 



/ 



dxe^ x (0\TfJx)j b »(0)\0 



lm[(0\Tj^q)j b »(-q)\0)}, 



(11) 



and expand the electromagnetic current (in momentum space) j^(±q) = (bj^b)(±q) which pro- 
duces/annihilates a bb pair with cm. energy y/q* in terms of 3 S± NRQCD currents up to dimension 
eight (i = 1,2,3) 



~ji(q) 



h (ftaix) (<?) 



6M t 2 



ji(-q) = h (x^iW) (-q) 



^a l (-^D) 2 X )(q) + ... , 

1,2 (xV i (-|B) 2 v;)(- g ) + ... ) 



6M t 2 



(12) 



where the constants 61 and b 2 are short-distance coefficients normalized to one at the Born level. 
Only the spatial components of the currents contribute contribute at the NNLO level. Inserting 



expansion fll2| ) back into Eq. (11) leads to the nonrelativistic expansion of the cross section at the 
NNLO level 



73thr 
NNLO 



(E) 



Ml 



Cl(Mhard)Mfac) Irn A\(E, /i so ft, Mfa 



3M, 4 



C^^hard; /if ac ) Im A 2 (E, ^ so ft, Mfac) 



+ 



(13) 



where 



Ax 
A 2 



U0\(ftax) (X^H D?^) +h.c.|0). 



(14) 
(15) 



The cross section is expanded in terms of a sum of absorptive parts of nonrelativistic current corre- 
lators, each of them multiplied by a short-distance coefficient. In fact, the right-hand side (RHS) of 
Eq. (FL3I) just represents an application of the factorization formalism proposed in O]. The second 



term on the RHS of Eq. (13) is suppressed by v 2 , i.e. of NNLO. This can be seen explicitly by using 
the equations of motion from from the NRQCD Lagrangian, which relates the correlator A2 directly 
to A\, 

A 2 = M b EA 2 . (16) 



Relation ( |16| ) has also been used to obtain the coefficient —4/3 in front of the second term on the 
RHS of Eq. (|b|). The nonrelativistic current correlators A\ and A 2 contain the resummation of the 
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singular terms mentioned in the previous paragraph. They incorporate all the long-distanceQ dynamics 
governed by soft scales like the relative three momentum ~ Mf,a s or the binding energy of the bb system 
~ The constants C\ and C2 (which are also normalized to one at the Born level), on the other 

hand, describe short-distance effects involving hard scales of the order of the bottom quark mass. 
They only represent a simple power series in a s and do not contain any resummations in a s . Because 
we consider the total bb cross section normalized to the point cross section, Eq. (|5|), C\ and C2 are 
independent of q 2 . In Eq. ( Jl3| ) we have also indicated the dependence of the correlators and the short- 
distance coefficients on the various renormalization scales: The factorization scale /if ac essentially 
represents the boundary between hard and soft momenta. The dependence on the factorization scale 
becomes explicit because of ultraviolet (UV) divergences contained in NRQCD. Because, as in any 
effective field theory, this boundary is not defined unambiguously, both the correlators and the short- 
distance coefficients in general depend on /Xf ac . The soft scale /U so f t and the hard scale /^hard) on the other 
hand, are inherent to the correlators and the short-distance constants, respectively, governing their 
perturbative expansion. If we would have all orders in a s and v at hand, the dependence of the cross 
section -R^nhd 011 variations of each the three scales would vanish exactly. (It is important that the soft 
and the hard scale, which both originate from the light degrees of freedom in the NRQCD Lagrangian, 
and the factorization scale are each considered as independent. They can each be defined in different 
regularization schemes. In this work we will use the MS scheme for the soft and the hard scale and a 
cutoff scheme for the factorization scale.) Unfortunately, we only perform the calculation up to NNLO 
in a s and v which leads to a residual dependence on the three scales /u.f ac , ^ so ft and /^hard- In particular, 
(as we will demonstrate in Section |||) the dependence on the soft scale /i so ft is quite strong, clearly 
because it governs the perturbative expansion of the correlators where convergence of the perturbation 
series can be expected to be worse than for the short-distance constants. It is therefore necessary to 
fix a certain window for each of the renormalization scales for which the perturbative series for the 
cross section shall be evaluated. At this point one can basically only rely on physical intuition, which 
tells that the renormalization scales should be of the same order as the physical scales governing the 
particular problem. This means that the soft scale should be the order of the relative momentum of 
the bb system^] ~ Mba s , and that the hard scale should be of order Mj, ~ 5 GeV. The factorization 
scale, on the other hand, should cover (at least partly) the soft and and hard regime. Because there 
is in our opinion no unique way to make this statement more quantitative, it is important to choose 
the corresponding windows "reasonably large" . In our case the choices are as follows: 

1.5 GeV < // soft < 3.5 GeV, 
2.5 GeV < ^hard < 10 GeV, 

2.5 GeV < /i fac < 10 GeV. (17) 

We will show in Sections |5| and ^ that the dependence of the theoretical moments P^ h on theses scales 
represents the dominant source of the uncertainties in the extraction of Thus, it is the choice 

4 In the context of this paper "long-distance" is not equivalent to "nonperturbative" . 

5 It is not clear at all whether there a not even smaller energy scales ~ Mb aj, k > 2, which might become relevant. 
However, those scales can only be produced by higher order effects like the hyperfine splitting, which should be irrelevant 
at least for the total cross section at NNLO. 

6 We will see later that at NNLO all interactions can be treated as instantaneous. As a consequence scales of the 
order of the binding energy ~ Mj, o? s can be ignored. 
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t t 

Figure 1: Graphical representation of the longitudinal and the transverse gluon exchange including the 
corresponding Feynman rules for the momentum exchange q = (q°,q). The exchange of a longitudinal 
gluon is instantaneous in time because its does not have an energy dependence. As a consequence the 
longitudinal exchange can be described by an instantaneous potential. The exchange of a transverse 
gluon, on the other hand, is retarded in time and, in general, cannot be described in terms of an 
instantaneous potential. 



given in Eqs. ( |17|) which determines the size of the uncertainties! 
Instantaneous Interactions and Retardation Effects 

To calculate the correlators A\ and A2 we use methods originally developed for QED bound state 
calculations in the framework of NRQED |?], 11, |l7|, |l8| and transfer them (with the appropriate 



modifications to account for the non-Abelian effects) to the problem of heavy quark-antiquark pro- 
duction in the kinematic regime close to the threshold. Because the Coulomb gauge is the standard 
gauge in which QED bound state calculations are carried out we also use the Coulomb gauge for the 
calculations in this work. The Coulomb gauge separates the gluon propagator into a longitudinal and 
a transverse piece (see Fig. |]). The longitudinal propagator does not have an energy dependence and 
therefore represents an instantaneous interaction. As a consequence, in configuration space representa- 
tion a longitudinal gluon exchange can be written as an instantaneous potential (which only depends 
on the spatial distance). Through the time derivative in the NRQCD Lagrangian the longitudinal 
gluon exchange leads to the Coulomb potential which is the dominant (LO) interaction between the 
bottom quarks in the nonrelativistic limit. Through the l/Mjf couplings of the bottom quarks to 
the chromo-electric field the longitudinal exchange also leads to the Darwin and spin-orbit potential, 
which contribute at the NNLO level. Because these potentials are instantaneous their treatment is 
straightforward in the framework of a two-body Schrodinger equation. 

For the transverse gluon the situation is more subtle. Because all couplings of the bottom 
quarks to the chromo-magnetic field are of order l/Mj the exchange of a transverse gluon between 
two bottom quark lines is a NNLO effect. However, in contrast to the Darwin and the spin-orbit 
interaction, the propagation of the transverse gluon energy has an energy dependence, i.e. it is an 
interaction with a temporal retardation. Physically this means that the transverse gluon can travel 
alongside the bb pair for some time period j|, . In this time period the bb pair is part of a higher 
order Fock 66-gluon state which, in principle, cannot be treated in terms of a two-body Schrodinger 
equation. Fortunately, in our case we can neglect the energy dependence of the transverse gluon 
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Figure 2: Graphical representation of the resummation of Coulomb ladder diagrams to all orders. 
The quark-antiquark propagation contains the nonrelativistic kinetic energy. The resummation is 
carried out explicitly by calculating the Green function the nonrelativistic Schrodinger equation with 
the Coulomb potential at the Born level, see Eq. (|27|). 




a b c 



Figure 3: Typical diagrams describing the exchange of a transverse gluons (in Coulomb gauge) in 
the background of the Coulomb exchange of longitudinal gluons. Longitudinal lines with a s ig n 
represent the summation of Coulomb ladder diagrams to all orders, see Fig. ^. 



propagator completely. This can be easily understood by considering a typical diagram describing 
the exchange of a transverse gluon between the bb pair in the background of a continuous Coulomb 
exchange of longitudinal gluons, see e.g. Fig. [3|a. If both ends of the transverse gluon end at bottom 
quarks the typical energy carried by the gluon can only be of order M^v 2 , the cm. kinetic energy of 
the bottom quarks. The typical three momentum of the gluon, on the other hand, can either be of 
order M&f, the relative momentum in the bb system, or also of order M^v 2 . If the three momentum 
is of order Mt,v 2 , the transverse gluon is essentially real and needs, in addition to the v 2 suppression 
coming from the couplings to the quarks, another phase space factor v to exist. Thus, the transverse 
gluons with this energy-momentum configuration lead to effects suppressed by v 3 , which is beyond the 
NNLO level. If the three momentum of the transverse gluon is of order Mj,v, on the other hand, it 
is far off-shell and we can neglect the small energy component in a first approximation. [It should be 
emphasized that this argument implies the hierarchy Mi,a s 3> M^a 2 , which is conceivable for the bb 
system where a s ~ O.3.] From that one can see that at NNLO the transverse gluon exchange can, like 
the longitudinal one, also be treated as an instantaneous interaction. This means that in Fig. ||a only 
those diagrams contribute at NNLO where the transverse line does not cross any longitudinal line. The 
differences between longitudinal and transverse gluons will only become manifest beyond the NNLO 
level. For the same reason any self-energy or crossed-ladder type diagram (see Figs. ||b,c for typical 
examples) can be safely neglected at the NNLO level. In fact, the situation is in complete analogy 
to the hydrogen atom or the positronium in QED, where it is well known that retardation effects 
lead to the "Lamb-shift" corrections which are suppressed by a 3 relative to the LO nonrelativistic 
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contributions. [Of course, the crossed exchange and self-energy type diagrams have to be taken into 
account in the two- loop calculation of the cross section in full QCD needed to determine the 0{a 2 s ) 
contributions to the short-distance coefficients. Those short-distance constants, however, describe 
effect from high momenta of order M b which are not contained in the correlators.] 

From the considerations above we can draw the following conclusions regarding the calculation 
of the correlators A\ and A2 at NNLO: 

1. We can treat the problem of bb production close to threshold as a pure two-body problem. This 
means that the NRQCD Lagrangian effectively reduces to a two-body Schrodinger equation from 
which the correlators can be determined. 

2. All interactions between the bottom and the antibottom quark can be written as time indepen- 
dent, instantaneous potentials, which means that only ladder-like diagrams have to taken into 
account. 



3. We can use the well known analytic solutions of the nonrelativistic Coulomb problem for positro- 



nium [20, 21, 22] and use Ray leigh- Schrodinger time- independent perturbation theory (TIPT) 
to determine the corrections caused by all higher order interactions and effects. 

However, there is one remark in order: although the effects of the transverse gluon exchange having 
a temporal retardation are formally beyond the NNLO level, this is not a proof that they are indeed 
smaller than the NNLO contributions calculated in this work. It is in fact rather likely that the 
retardation effects cannot be calculated at all using perturbative methods because the characteristic 
scale of the coupling governing the emission, absorption or interaction of a gluon which has energy 
and momentum of order M b a 2 s would be of the order of 0.5 — 1 GeV. This is already quite close to the 
typical hadronization scale A QCD . From this point of view it seems that the NNLO analysis presented 
here cannot be improved any more, at least not with perturbative methods. This problem might even 
cast doubts on the reliability of the NNLO corrections themselves and underlines the necessity that 
the preferred ranges for the renormalization scales, Eqs. (0), are chosen sufficiently large. We will 
ignore further implications of this problem for the calculations and analyses carried out in this work. 
(See also the paragraph on nonperturbative effects.) 



Instantaneous Potentials 

At the Born level all potentials relevant for the nonrelativistic cross section at NNLO can be obtained 
directly from the NRQCD Lagrangian considering (color singlet) bb —* bb single gluon t-channel 
exchange scattering diagrams. In configuration space representation the Born level potentials read 
(r = \r\, C F = 4/3, C A = 3, T = 1/2, a s = a s (// soft )) 

Cf a s 



r 

C F a s 7r 

~w 

3 Cf a 



1 + -S b Si 



Cp a s 
2M b 2 r 



v- 



+ -^r (fV)V 



M^r 3 



1 1 



(18) 



where S b and S b are the bottom and antibottom quark spin operators and L is the angular momen- 
tum operator, vj ^ is the well known Coulomb potential. It constitutes the LO interaction and will 
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(together with the nonrelativistic kinetic energy) be taken into account exactly. It arises from the 
exchange of a longitudinal gluon through the time derivative coupling of the bottom quarks to the 
gluon field. Vbf represents the Breit-Fermi potential which is known from higher order positronium 
calculations. It describes the Darwin and spin-orbit interactions which are mediated by the longitudi- 
nal gluons and also the so called hyperfine or tensor interactions which are mediated by the transverse 
gluons in the instantaneous approximation. Due to the 1/-M? suppression V BF already leads to NNLO 
effects in the cross section and will be taken into account as a perturbation. For the same reason 
only the radiative corrections to the Coulomb exchange of longitudinal gluons have to be taken into 
account. We want to emphasize that these radiative corrections are caused by the massless degrees of 
freedom in the NRQCD Lagrangian. Because in the corresponding loops transverse gluon lines can end 
at other massless lines the considerations given in the preceding paragraph cannot be applied in this 
case. Thus, in general, transverse gluons (or massless quarks) in all energy-momentum configurations 
have to taken into account to calculate the radiative corrections properly. The calculation of these 
radiative corrections can be found in existing literature and we therefore just present the results. 

At the one-loop level (and using the MS scheme for the strong coupling) the corrections read 
(7e = 0.57721566 . . . being the Euler-Mascheroni constant) 



where 
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m 
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2/3 ln(/2r) + a\ 



4 

-7 

3 
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A = e 7E ^ soft , 



(19) 



(20) 



and 



C A C F a 2 s 
2M b r 2 



(21) 



y c (1) 



represents the one- loop corrections to the Coulomb potential oc 1/r and leads to NLO contri- 
butions in the cross section, vi 1 ^ has been calculated by Fischler [23] and Billoire V^ A) called 
non-Abelian potential for the rest of this work, arises from a non-analytic behavior of the vertex 
diagram depicted in Fig. |] oc (k 2 /M 2 ) 1 / 2 where k is the three momentum exchanged between the 
bottom and the antibottom quark. Because the non-analytic term causes a behavior oc l/\k\ for the 
non-Abelian potential in momentum space representation, V NA is proportional to 1/r 2 . We would like 
to point out that in Coulomb gauge such a non-analytic behavior does not exist for Abelian diagrams. 
We refer the reader to [^5[ for older publications, where the non-Abelian potential has been deter- 
mined. Due to the a s jM\, factor V^ A is a NNLO interaction and no further corrections to it have to 
taken into account. 

At the two-loop level only the corrections to the Coulomb potential have to be considered. 
They have been calculated recently by Peter |2(| and read (in the MS scheme) 

2 



K (2) (r 



4 1n 2 (/ir) + ^- j + 2(2/3 ai + /3i) ln(/i r) + a 2 



(22) 
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Figure 4: Vertex diagram in Coulomb gauge responsible for the potential non-Abelian potential V NA . 
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Figure 5: Symbols describing the interactions potentials Vc°\ Vc , V c , V BF and V NA and the kinetic 
energy correction 5H kin = -V 4 /4M 6 3 . 



where 



Pi = Y c A-f CATni ~ ACFTni > 

( 4343 9 vr 4 22 \ o / 1798 56 \ „ m 

5-16Cs)cFT ni + (frm) 2 . (23) 

For later reference we assign the symbols in Fig. [5| to the potentials given above. We also would like 
to note that we do not have to consider any annihilation effects. The leading annihilation diagram is 
depicted in Fig pi Because the annihilation process takes place at short distances, it produces local 
four-fermion operators in the NRQCD Lagrangian, which can be written as instantaneous potentials. 
The dominant annihilation potential which comes from the three gluon annihilation diagram has the 
form Va nn (r) oc (a^/M^)5^ (r) and would lead to effects suppressed by v A in the cross section. 

Recipe for the Calculation of Large n Moments at NNLO 

Based on the issues discussed above the calculation of the NNLO nonrelativistic cross section -R*nlo 
and the theoretical moments P^ h in terms of the correlators A\ and A2 and the short-distance coeffi- 
cients C\j2 proceeds in the following three basic steps: 
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<mns> 



<nrmrp> 
nmnnp 

Figure 6: The dominant annihilation diagram relevant for bb — > bb scattering for a bottom-antibottom 
quark pair in a color singlet J PC = 1 , 3 S\ configuration. Its dominant contribution leads to a 
potential Knn(^) oc a? s /M 2 5(f) and to contributions in the cross section and the moments beyond the 
NNLO level. 



Step 1: Solution of the Schrddinger equation. - The Green function of the NNLO Schrodinger 
equation is calculated incorporating the potentials displayed above and including the 
NNLO corrections to the kinetic energy. The correlators A\ and A2 are directly related 
to the zero-distance Green function of the Schrodinger equation. 

Step 2: Matching calculation. - The short-distance constant C\ is determined at 0(a 2 s ) by match- 
ing expression (|l3|) directly to the cross section calculated in full QCD at the two-loop level 
and including terms up to NNLO in an expansion in v in the (formal) limit a 8 -C ti < 1. 

Step 3: Dispersion Integration. - The integration (Q) is carried out. 

For the rest of this section we briefly explain the strategies and basic procedures for the steps 1 and 2. 
The explicit calculations for steps 1-3 are presented in detail in Sections 3.1, 3.2 and|3.3|, respectively. 



Solution of the Schrddinger equation: The nonrelativistic correlators .4.1 and A2 are calculated by 



determining the Green function of the Schrodinger equation (E 
V 4 



2M b ) 



+ 



V W(jf) + V P(f) + V^(f) + V BF (f) + V NA (f) 

= 5 {3) (f-f 



E^j G(r,r*,E) 



(24) 



where Vbf is evaluated for the 3 S± configuration only. The relation between the correlator .Ai and 
Green function reads 



Ai 



6 N r . 



lim G(f,f,E) 



(25) 



Eq. (|25| ) can be quickly derived from the facts that G(f,f f ,E) describes the propagation of a bottom- 
antibottom pair which is produced and annihilated at relative distances \r\ and \ / r i \, respectively, 
and that the bottom-antibottom quark pair is produced and annihilated through the electromagnetic 
current at zero distances. Therefore A\ must be proportional to limij?i ifi— >o G(f,r l ,E). The correct 
proportionality constant can then be determined by considering production of a free (i.e. a s = 0) 
bottom-antibottom pair in the nonrelativistic limit. (In this case the Born cross section in full QCD can 
be easily compared to the imaginary part of the Green function of the free nonrelativistic Schrodinger 
equation.) The correlator A2 is determined from A\ via relation (|l~6|). We would like to emphasize 
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that the zero-distance Green function on the RHS of Eqs. (Ejj) contains UV divergences which have to 
regularized. In the actual calculations carried out in Section |3.l| we impose the explicit short-distance 
cutoff jUfac- As mentioned before, this is the reason why the correlators and the short-distance constants 



depend explicitly on the (factorization) scale /if ac . In this work we solve equation (24) perturbatively 



by starting from well known Green function Gc of the nonrelativistic Coulomb problem |29, 21 



V 2 

TT ~ V c°\r) ~ E ) G c (?y,E) = 6®(?- j*) (26) 



M b 

and by incorporating all the higher order terms using TIPT. 

Matching calculation: After the nonrelativistic correlators Ai and A2 are calculated the determina- 
tion of C\ is achieved by considering the (formal) limit a s < » < 1. In this limit fixed multi-loop 
perturbation theory (i.e. an expansion in a s ) as well as the nonrelativistic approximation (i.e. a sub- 
sequent expansion in v) are feasible. This means that multi-loop QCD (with an expansion in v after 
the loop integrations have been carried out) and multi-loop NRQCD must give the same results. In 
our case we use this fact to determine the constant C\ up to terms of order a 2 . For that we expand 
the NNLO NRQCD expression for the cross section ( |l~3| ) for small a s up to terms of order a 2 and 
demand equality (i.e. match) to the total cross section obtained at the two-loop level in full QCD 
keeping terms up to NNLO in an expansion in v. Because NRQCD is an effective field theory of QCD 
(i.e. it has the same infrared behavior as full QCD) for the limit v -C 1, Ci contains only constant 
coefficients (modulo logarithms of the ratios Mb/fi^c arid Mft/jUhard)- All the singular terms oc 1/vAivv 
are incorporated in the correlators Ai and A2 



Comment on Nonperturbative Effects 

To conclude this section we would like to mention that nowhere in this work nonperturbative ef- 
fects in terms of phenomenological constants like the gluon condensate ( | G^y G^ v | ) [lf|] are taken 
into account. In ^] it has been shown that the contribution of the most important condensate 
( \Gfj, v G^ v \ 0) is at the per-mill level in the moments P™ 1 for 4 < n < 10. As we show in Sec- 
tion ^[ this effect is completely negligible compared to the theoretical uncertainties coming from the 
large renormalization scale dependences of the NNLO moments P^ h . The condensates are therefore 
irrelevant from the purely practical point of view. 

Nevertheless, we even think that the inclusion of the condensates for the moments at the NNLO 
level would be conceptually unjustified. For the gluon condensate this can be seen from the fact that it 
provides a phenomenological parameterization of the average long-wavelength vacuum fluctuations of 
the gluon field involving scales smaller than the relative three momentum of the bb system Thus, 
for the theoretical moments P„ (4 < n < 10) (and also for heavy enough quarkonia in general) the 
condensates describe retardation-like effects ||. As explained before, we neglect retardation effects 
because they formally contribute beyond the NNLO level. We conclude that taking into account the 
condensates would only be sensible in a complete NNNLO analysis. In this respect the condensate 
contributions might provide some estimates for the size of some NNNLO effects. However, if the small 
size of the condensate effects in the moments P^ 1 is compared to the large perturbative uncertainties 
contained of the NNLO theoretical moments, it is seems rather doubtful whether the condensates 
represent the dominant contributions at the NNNLO level. 
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3 Calculation of the Moments 



In this section the determination of the theoretical moments is presented in detail. Because all 
conceptual issues have been discussed is Section [j] we concentrate only on the technical aspects. The 



task is split into three parts which are described in the following three subsections. In Section 3.1 



the nonrelativistic correlators A\ and A2 are calculated and Section 3.2 describes the calculation of 
the short-distance constant C\. In Section [3l^ the dispersion integration (B) is carried and the final 
formulae for the theoretical moments are presented. 

3.1 Calculation of the Nonrelativistic Correlators 

To calculate the nonrelativistic correlators A\ and A2 the Green function G of the Schrodinger equa- 
tion ( p4| ) has to be determined. As explained before, we start with the Green function of the 
nonrelativistic Schrodinger equation (Efj), called "Coulomb Green function" from now on, and de- 
termine the effects all the higher order contributions through TIPT. The most general form of the 
Coulomb Green function reads (r = \r\, r' = jr"!) 



G^(fy,E) 

d 2 

X dtds 



M h 



4 7iT(l + ip)r(l-i/j) 



1 00 
dt I ds 



s(l-t) t(s-l) 



-ip 







1 



ts 



sr — tr 



exp < ip (^\r \ (1 — t) + \f\ (s — 1) + \sr — t'i 



r < r . 



where 



p = M b v = J M b (E + ie), 



P 



Cf a s 
2v 



(27) 



(28) 



and r is the gamma function. The case r < r' is obtained by interchanging r and r'. G^\f,r',E) 
represents the analytical expression for the sum of ladder diagrams depicted in Fig. ^. We refer the 
reader interested in the derivation of G^ to the classical papers 21, 22]. The analytic form of 
the Coulomb Green function shown in Eq. ( p7| ) has been taken from Ref. |2C]. Fortunately we do not 
need the Coulomb Green function in its most general form but only its S-wave component 



G*°) ' S (r, r' , E) = ^- I dQ (r , f , E) 

47T J 

2iM b p 



4irT(l + ip)T(l-ip) 



1 00 
dt / ds 



ip 



s(l-t) t(s-l) 



-ip 







1 



exp|ip r'(l-2t)+r(2s-l) }, 



r < r . 



(29) 



The case r < r' is again obtained by interchanging r and r' . For r' = the form of the Coulomb 
Green function is particularly simple 



G[ 0) (0,r,£) = G^(0,r,E) 
.M b p 



M b p 
2vr 



dte 



2iprt 



1 + t 



1 P 



2n 



e ipr T(l-ip)U(l-ip,2,-2ipr) 



15 




4irr 



T(l-ip)W ipt i(-2ipr). 



(30) 



where U(a,b,z) is a confluent hypergeometric function and W K u(z) one of the Whittaker's func- 
p8| . It is an important fact that G^ (0, r, E) diverges for the limit r — > because it contains 



tions ||27], p 

power (oc 1/r) and logarithmic (oc lnr) divergences As explained in Section |l| these ultraviolet 
(UV) divergences are regularized by imposing the small distance cutoff /if ac . The regularized form of 
lim r ^o G c (0,f,E) reads 



/If 2 f 

GP' re9 (0,0,E) = -JL\iv-C F a s 

4 7T ^ 



ln(- 



.M b v, 

i 

Piac ' 



+ Te + * 1 



2v 



(31) 



where the superscript "reg" indicates the cutoff regularization and ^(z) = dlnT(z) /dz is the digamma 
function. For the regularization we use the convention where all power divergences oc ;Uf ac are freely 
dropped and only logarithmic divergences oc ln^fac/Mf,) are kept. Further, we define /if ac such that 
in the expression between the brackets all constants except the Euler-Mascheroni constant 7 E are 
absorbed. The same convention is also employed for the calculation of the higher order corrections 
to the Coulomb Green function which are discussed below. The results for any other regularization 
scheme which suppressed power divergences (like the MS scheme) can be obtained by redefinition of 
the factorization scale. Our apparently sloppy realization of the regularization procedure is possible 
because in Section we will match the expression for the NNLO cross section in NRQCD directly 
to the corresponding two-loop expression in full QCD. As a consequence additional constant terms 
in the brackets on the RHS of Eq. ([}(]) do not affect the final result for the cross section at NNLO 
in NRQCD because they merely represent contributions which can be anyway freely shifted between 
the nonrelativistic correlators and the short-distance coefficients. The resulting (small) ambiguity 
will be accounted for during the fitting procedure by varying the factorization scale /if ac - However, 
the reader should note that even with a more stringent regularization scheme like the MS scheme 
this ambiguity cannot be avoided because the factorization //f ac is by construction not fixed to any 
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value regardless which regularization scheme is used. 7 For later reference we call G^' reg (0,0,E) 
"zero-distance Coulomb Green function". A graphical representation of G^' re9 (0, 0, E) in terms of 
NRQCD Feynman diagrams is displayed in Fig. ^a. For convenience we suppress the superscript "reg" 
from now in this work. 

The Coulomb Green function contains bb bound state poles at the energies \fs n = 2M& — 
CpCL^Mb/An 2 (n = 1, 2, . . . oo). These poles come from the digamma function in Eq. ( |3l| ) and corre- 



spond to the nonrelativistic positronium state poles known from QED [29|. They are located entirely 
below the threshold point y^thr = This can be seen explicitly from the cross section in the 

nonrelativistic limit, 



^Im 

Ml 



6TrN c Q 2 b 
o — - 1m 

Ml 



G<°>(0,0,£) 



24vr 2 iV c <3 



M b 



oo 

E 

n=l 



+ ®{E)-N c Ql 



Cp a s 7T 



1 — exp(- 



(32) 



where ^^(O)! 2 = (MbCFa s ) 3 /8nn 3 is the modulus squared of the LO nonrelativistic bound state wave 
functions for the radial quantum number n. The continuum contribution on the RHS of Eq. ( |32|) is 
sometimes called "Sommerfeld factor" of "Fermi factor" in the literature. The resonance contributions 
are described by the first term in the second line of Eq. (|32|). The corrections to the zero-distance 
Coulomb Green function calculated below lead to higher order contributions to the bound state energy 
levels, the residues at the bound state poles and the continuum. We would like to stress that all these 
contributions must be included in the dispersion integration (Q) to arrive at reliable results for the 
theoretical moments P„ . Nevertheless, it is worth to note that the resonances are not necessarily 
equivalent to the actual T resonances |J . In particular for large radial excitations a direct comparison 
would be more than suspicious. In the context of the calculation of the moments they have to be 
included for mathematical rather than physical reasons. (See also the discussion in Section H[) 

Let us now come to the determination of the corrections to the zero-distance Coulomb Green 



function coming from the remaining terms in the Schrodinger equation (24). At NLO only the one- 



loop contributions to the Coulomb potential, vj 1 ^ (see Eq. (|19|)), have to considered. Using first order 
TIPT in configuration space representation the NLO corrections to Gc (0,0, E) read 



G«(0,0,£) 



d 3 f G^ (0, r, E) V c m (r) G^ (r, 0, E) . 



(33) 



Expression (|33|) is displayed graphically in Fig. [?]b. Further evaluation of the integration on the 
RHS of Eq d33| ) is possible but not presented here, because it is already in a form suitable for the 
dispersion integration (0) (see Section |3.3[ ). At NNLO several contributions have to be considered. 
The corrections from the two-loop contributions to the Coulomb potential, Vc" 2 ^ (see Eq. (p^)), are 
calculated in analogy to the NLO contributions using first order TIPT (Fig. ^c) 



Gi 2) (0,0,2*0 



2 loo 



d 3 f G^ (0, r, E) y c (2) (r) G^ (r, 0, E) . 



(34) 



We also have to take into account the one-loop Coulomb potential (Fig. |7|(d)) in second order TIPT, 



G^(0,0,E) 



1 loo 



7 In fact, using the MS scheme is a quite tricky (but not impossible) task if one wants to avoid solving the Schrodinger 
equation in D dimensions. 
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d 3 ri / d 3 f 2 (0, n , £0 ft) Gi 0),S (n, r 2 , £) ^i 1 ) ft) G® (r 2 , 0, -E) . 



(35) 



Because the Coulomb potential is angular independent, only the S-wave components of the Coulomb 
Green function in the center of expression ( |35|) are needed. Finally, we have to determine the NNLO 
contributions to the zero-distance Green function coming from the kinetic energy, 5H kin = — V 4 /4M^, 
the Breit-Fermi potential V BF and the non-Abelian potential V NA (see Figs, flje and f). These correc- 

( 2) 

tions are symbolized by [Gc (0, 0, -E)] km+BF+NA in the following. A method to determine them has 
been presented in an earlier publication |7|, |1Q]. Some details about this method are presented in 
Appendix The final result for [Gi 2) (0, 0, £)] kin+BF+NA reads 



G^(0,0,E) + \g {2) (0,0,E) 
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2v 



Cf Qj 
2v 



(36) 



(2) 

Because [Gc (0, 0, £')] km + BF + NA contains also kinematic corrections to the zero-distance Coulomb Green 
function, we found it convenient to add the zero-distance Coulomb Green function (|3~I|). The first 
term on the RHS of Eq. ( |36|) represents the zero-distance Coulomb Green function including the 
NNLO kinematic corrections and the second term the remaining corrections. It is an interesting fact 
that these remaining corrections can be written as the squared of the zero-distance Coulomb Green 
function. This is a consequence of the (renormalization group) invariance of the total cross section (13) 
under variations of the factorization scale /i S oft- (See also the comment after Eq. (|64|).) Collecting all 
contributions the complete expression for the nonrelativistic correlator Ai at NNLO reads 



= 6N C { G(°)(0,0,£)+G«(0,0,£) 
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4 2) (0,0,£) 



1 loo 
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G( 2) (0,0,£) 
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(37) 



The calculation of the correlator A2 , on the other hand is trivial using the equation of motion for the 
Green function, see Eq. (|i~6|). Because A2 is multiplied by an explicit factor v 2 , Eq. (|i~6|), its form is 
particularly simple, 
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3M£ 
2vr 



1 v 



Cf a s 
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+ 7E + ^(l 



2v 



(38) 



3.2 Determination of the Short-distance Coefficients 

The short-distance coefficient C\ and C 2 are determined by matching the NNLO cross section (|b|) in 
NRQCD to the same cross section calculated in full QCD (in the limit a s <C v <C 1) at the two-loop 
level and including terms in the velocity expansion up to NNLO. It is convenient to parameterize the 
higher order contributions to C\ in the form (ah = a s (/^hard)) 



Ci(M b ,^ hard ,/i fac ) = 1 + (^) 4^ + (-^) cfVhard,Alfac) 



+ ... . 



(39) 
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Due to renormalization group invariance only the 0(a 2 ) coefficient of C\ depends on the hard scale 
^hard- We have already anticipated that the 0{a s ) coefficient does not depend on the factorization 
scale //f ac . For C2, on the other hand, no higher order contributions are needed because the correlator 
A 2 is already of NNLO, 

C 2 = 1. (40) 

The expansion of the NNLO cross section in NRQCD, -R^nlo' ^q. (13), keeping terms up to order a 2 , 
reads 



Dthr C*sjCl tit p. 
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192 



+ 0(«s 



(41) 



where we have set /i S oft = Mhard because in the limit a s C » < 1 a distinction between soft and hard 
scale is irrelevant. (We want to emphasize that the choice /^ so ft = A*hard implicitly means that the hard 
scale is, like the soft scale, defined in the MS scheme.) The corresponding expression for the two-loop 
cross section calculated in full QCD readsPI 
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(42) 



+ <V) , 
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(43) 



The Born and one-loop contributions in Eq. (|42| ) are standard [.BC], 31 1. The two-loop contributions 
are presented with the various combinations of the SU(3) group theoretical factors Cf = 4/3, Ca = 3 
and T = 1/2. The terms proportional to Cp come from the QED-like, Abelian exchange of two gluons 



and have been calculated analytically in |32|]. The result has been confirmed numerically in [33| and 
analytically in |34|, ^5|. The corresponding Feynman diagrams (in the covariant gauge) are displayed 
in Fig. ||](a), (b), (c) and (d). The CaCf terms correspond to the non-Abelian exchange of two gluons, 

8 The two-loop contributions from secondary radiation of a bb pair off a light quark-antiquark pair are kinematically 
suppressed and do not contribute at NNLO in the velocity expansion. 
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Figure 8: QCD Feynman diagrams relevant for the calculation of the cross section at the two-loop 
level. The calculation of these diagrams is needed for the matching calculation which leads to the 
determination of the short-distance coefficient C\. Feynman diagrams needed for the wave function 
renormalization are not displayed. 



i.e. involving the triple gluon vertex, ghost fields and topologies with crossed gluon lines (Figs. ||](b), 
(c), (d), (e), (f), (g)). These contributions have been determined in |34|, |35|. The CfThi contributions 
are from diagrams with a vacuum polarization of massless quarks (Fig. ||(h)) and have been calculated 
in ||36|| . The contributions proportional to CjrT, finally, correspond to the diagram where the vacuum 
polarization is from the bottom quarks (Fig. ||(g)) and have been calculated in |37] , |3(|. The virtual 
top quark contributions are suppressed by a factor (M^/Mt) 2 ~ 0.001) and are neglected. 

The constants and defined in Eq. (|39| ) can now be easily determined by demanding 
equality of expressions (fH] ) and (f42|). This constitutes the "direct matching" procedure [0, [ll| and 
leads to 
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(44) 
(45) 



The constant Cj 1 is the 0(a s ) short-distance contributions which is well known from the single photon 
annihilation contributions to the positronium hyperfine splitting [|38| and from corrections to electro- 
magnetic quarkonium decays |}9| . We want to mention again that //hard and ^f ac are independent and 
defined in different regularization schemes. 

To conclude this subsection we would like to point out that the short-distance coefficients C\ 
and C2 determined above are not sufficient to determine the vacuum polarization function (Eq. (|])) 
in the threshold regime at NNLO, because they have been determined via matching at the level of 
the cross section only, i.e. at the level of the imaginary part of the vacuum polarization function. 
The expressions for the correlators still contain overall UV divergences oc ln(M(,//if ac ) in their real 
parts [19, 29], see e.g. Eq. (31). For the large n moments calculated in this work these ambiguities 
are irrelevant because the divergent contributions in the real parts do not contribute to the large n 
moments. The relation between the nonrelativistic correlators and the vacuum polarization function 
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at NNLO in the threshold regime, including the proper short-distance contributions for the real part, 
has the form 
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The constants /ii and /12 can be determined via (direct) matching to the one and two-loop vacuum 
polarization function in full QCD at threshold, i.e. for q 2 — > AM 2 . This work has been carried out in a 
previous publication |l6| and leads to h\ = and /12 = 4^7 (3 — ^"Cs) + 55 ~~ f m 2- For the complete 
expression of the vacuum polarization function in the threshold regime at NNLO in the nonrelativistic 
expansion also the 0(a 2 ) and 0(a 3 ) short-distance contributions would have to be calculated. This 
would require the calculation of the three- and four loop the vacuum polarization functions in full 
QCD in the threshold regime. This task has not been accomplished yet and remains to be done.^] 

3.3 The Dispersion Integration 

After the nonrelativistic correlators A\ and A% and the short-distance constants C\ and C2 are calcu- 
lated we are now ready to carry out the dispersion integration (|j). This task is quite cumbersome if 
the complete covariant form of the integration measure ds/s n+1 is used. Fortunately the integration 
can be simplified because we are only interested in NNLO accuracy in the nonrelativistic expansion in 
v = (E/Mb) 1 / 2 . Changing the integration variable to the energy E = y/s — 2M b and expanding up to 
NNLO in v, where combination (E/M b )n is considered of order one, the resulting integration measure 
reads 

ds 1 dE f . , ,/ E 

exp <^ - (2 n + 1) In 1 + 



(4M 2 ) n M b y \ y ' ; V 1 2M b 
E^M b 1 dE [ E \ ( E E 2 E 2 E 3 E 4 2 



(4 M 2 )« M b 6XP 1 M b n )V 2 M b + 4 M 2 " + ° { M 2 ' M 3 M 4 ^ 



The dispersion integration for the theoretical moments at NNLO then takes the form 



(47) 



p " = jdfr I i exp { " iH ~ tm, + £i;- n > R ™- (£) - <48> 



u.l 



where .Ebind is the (negative) binding energy of the lowest lying resonance. We would like to point out 
that expansion ( f47| ) leads to an asymptotic series, which means that including more an more terms 
in the expansion can improve the approximation only up to a certain point beyond which the series 

9 In Refs. [^(J numerical approximations for the three loop vacuum polarization valid for all energies has been obtained 
based on the Pade method. Unfortunately numerical approximations are of little use for a precise extraction the 0(o? a ) 
short-distance constants due to the presence of singular terms cx In v and In 2 v in the real part of the three loop vacuum 
polarization function close to the threshold. 
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Figure 9: Path of integration to calculate expression (48) for the theoretical moments P*' 1 . The 
dashed line closes the contour at infinity and does not contribute to the integration. The constant 7 is 
chosen large enough to be safely away from the bound state poles which are indicated by the gray dots 
on the negative energy axis. The thick gray line on the positive energy axis represents the continuum. 



starts diverging. We have checked that for all values of n employed in this work the expansion is still 
well inside the converging regime. It should also be noted that for increasing values of n the expansion 
provides better and better approximations only as long as the condition (E^ in( ^/M b )n < 1 is satisfied. 
In our case, where the bb system is treated as Coulombic, i.e. -Ebind = M b Cpa 2 /4: + .... this condition 
is always satisfied. (See also discussion at the end of Section |4[) Integration (|48"D is carried out most 
efficiently by deforming the path of integration into the negative complex energy plane as shown in 
Fig. ||[ Because the (dashed) line which closes the contour at infinity does not contribute and because 
we take 7 large enough to be safely away from the bound state poles (7 S> £?bind)j we can rewrite 
expression 



as 



P 



tli 



-2iQ 2 b 7T 

(4M b 2 )"+ 1 



-7+100 

2 - ' dE 



-7—400 

7+ioo 



1 



dE 



cxp 



E 



M~ h n ) V 2M b ' 4M 2 



E E 2 



E W E t E 2 n 
(4M b 2 ) n + 1 2vri J M b eW {W b n j \ + 2M~ b + AM 2 " 

7— too 



CiAx(E) - —C 2 A2(E) 



3M 2 



C 2 A 2 (-E) 



(49) 



where in second line the change of variables E — * —E has been performed. The reader should note 
that for the integration in the negative complex energy plane also the real part of the correlators *4i 
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and A2 is needed. The expression in the second line of Eq. (|49| ) offers three advantages which make 
it much easier to calculate than expression ([48|): 

1. Because the integration path is far away from bound state energies, the integrand can be ex- 
panded in a s . This avoids that we have to integrate over complicated special function like the 
digamma function ^. 

2. We do not have to integrate separately over the resonances and the continuum. Both contribu- 
tions are in a convenient way calculated at the same time. 

3. The expression in the second line of Eq. (|4^) is nothing else than an inverse Laplace transform 
for which a vast number of tables exist in literature (see e.g. [pqj ). 

We want to stress that the advantages described above are merely technical in nature and just simplify 
the calculation. The results of the integration are not affected. 

The final result for the theoretical moments including all contributions up to NNLO in the 
nonrelativistic expansion can be cast into the form 

3N c Q 2 b 



P. 



th 



7T 



Cl^hardj Mfac) Qn,l (^soft , Wac) + C*2 Qn,2 



(50) 



" 4(4M & 2 )™n 3 / 2 

where g n i comes from the integration of the correlator Ai (including LO, NLO and NNLO contribu- 
tions in the nonrelativistic expansion) and g nt 2 originate from the integration of A2 which is of NNLO 
only. To illustrate the technical aspects of the integration ( [49|) let first present some of the details 
of the calculation of the LO contribution to Q n \. The LO contributions to g n \ originates from the 
zero-distance Coulomb Green function in Eq. (|3~l|). The corresponding integration takes the form 
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C F a s In v + C F a s ^ ( p 

p=2 



2v 



where 




(51) 



(52) 



and Q v is the Riemann zeta function for the argument p. Because \Cfcl s /2v\ <C 1 along the integration 
path we have expanded the digamma function in G^(0, 0, — E) for small a s . The resulting expression 
is now immediately ready for the application of inverse Laplace transforms. Here, we only need the 
relations 



u-l 



7+ioo 
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7— ioo 
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(53) 
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The result for [g n i] LO reads 



i LO 



Qn.l 



p=2 



where 



i + 2 v / ^ + 4 v / ^£ </> p -f 



Expression ( |54|) can be rewritten in the form 

2vr 2 



Qn,l 



l + 2^(f> + 



p=i 



P 



1 + erf 



(54) 



(55) 



(56) 



where erf is the error function defined as erf (z) = ^= /q exp(— t 2 )dt. Expression d56| ) agrees with the 
result obtained by Voloshin ||. The infinite series defined in Eq. ( |54| ) is absolute convergent with 
an infinite radius of convergence. For the values of n employed in this work (4 < n < 10), however, 
convergence is somewhat slow and a large number of terms have to be taken into account. This fact is 





0.5 


0.6 


0.7 


0.8 


0.9 


1.0 


[£n,l] L ° 


6.38 


9.44 


14.07 


21.16 


32.10 


49.12 


first three terms in Eq. (p3|) 


4.42 


5.50 


6.71 


8.05 


9.52 


11.12 



Table 1: Comparison of the series for [^ n ,i] LO with the sum of Born, one- and two-loop contributions 
in the series on the RHS of Eq. fl54]) for the values of <fi employed in this work. 



illustrated in Tab. [I] where the sum of the first three terms (corresponding to Born, one- and two-loop 
contributions) in the series ([)4]) is compared to the total sum for values of (j> between 0.5 and 1.0, 
which represent the range of <f> values used in this work. Tab. |l] shows that the resummation of higher 
orders in a s is essential to arrive at sensible results in particular for larger values of n. This feature 
remains true for all contributions to g n i and g n ^ and shows that a naive fixed order (multi-loop) 
calculation for the moments is unreliable for large values of n. 

Along the lines of the calculation of [^ n ,i] LO it is now straightforward to determine g n 2 and the 



NLO and NNLO contributions to g 



corrections to the Coulomb potential, V^ 1 and vj 2 ^ 



nt l. The contributions to g n ^\ coming from the one- and two-loop 
have the form 
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The coefficients of the beta function, /3q,i and the constants ai^ are given in Eqs. ( pOj ) and (23). 
The function VP' is the derivative of the digamma function and 0F2 is a generalized hypergeometric 



25 



function [28]. The constants iv®' 1 ' 2 and w^ 1,2 are given in Appendix [C]. For the calculation of expres- 
sion (||^) the table of inverse Laplace transforms given in Appendix |b] has been used extensively. The 



term proportional to 5\ in Eq. (57) contains the NLO contributions coming from V c (1) and the NNLO 

/o\ 

contributions coming from the terms oc 1/r and oc ln(/x so f t r)/r in Vc in first order TIPT. The term 

proportional to 82 contains the remaining NNLO corrections coming from the term oc ln 2 (// so f t el r)Jr 
(2) 

m Vr>. The expression proportional to 03, finally, arises from the second order interaction in TIPT of 
Vc (1) - The NNLO contributions to g n 1 originating from the kinetic energy corrections, the Breit-Fermi 
potential, the non-Abelian potential (see Eq. ( |36|) for the corresponding corrections to the zero-distance 
Green function) and the kinematic correction factor (1 + E/2M{, + (E 2 /4M 2 )n) from Eq. ([48] ) read 
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(64) 



where, for convenience, also the LO result from Eq. (^3|) has been added. From the last line of expres- 
sion (^) one can easily determine a renormalization group equation (with respect to the factorization 
scale //f ac ) for g n> i and g n ,2, which would allow for a resummation of the corrections coming from the 
kinetic energy corrections, the Breit-Fermi potential and the non-Abelian potential to all orders in 
TIPT. Although it is quite tempting to carry out this resummation, we refrain from doing so because a 
resummation of those corrections would not account for the retardation effects mentioned in Section [| 
The complete expression for g n 1 has the form 
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(65) 



Finally, the result for g n ^ coming from the integration of A2, Eq. (pq), reads 
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(66) 



From from expression ( |50| ) for the theoretical moments at NLO one can easily recover the 
moments at NNLO by setting 
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Moment 




a s (M z ) 




4.6 


4.8 


5.0 


5.2 


0.10 


0.11 


0.12 


0.13 


Pl h /\W- xx GeV~ 8 l 


0.51 


0.37 


0.27 


0.20 


0.19 


0.27 


0.41 


0.74 


P£ h /[10- xx GeV" 12 ] 

Oil J 


0.67 


0.41 


0.25 


0.16 


0.17 


0.26 


0.46 


0.97 


P^ h /[10- xx GeV~ 16 ] 


0.95 


0.49 


0.26 


0.14 


0.18 


0.29 


0.57 


1.37 


pth /[IQ-xx GeV" 20 ] 


1.42 


0.61 


0.27 


0.13 


0.19 


0.34 


0.73 


1.99 


P 2 ^/[10- ra GeV- 40 ] 


12.96 


2.37 


0.47 


0.10 


0.42 


1.00 


3.07 


13.93 




a.{M z ) = 0.118 


M b = 4.8 GeV 




/"soft = 2.5 GeV , ^hard = Mfac = 5 GeV 



Table 2: The theoretical moments for n = 4,6,8, 10,20 and fixed // so ft = 2.5 GeV and /^hard = 
Wac = 5 GeV for various values of and a s (M z ). The two- loop running for the strong coupling has 
been employed. 

62 = 5 3 = , 

Mi = Msoft exp ^ ^ , (67) 

and by ignoring the corrections [£n,i]kin+BF+NA- The resulting expression for the NLO moments is 
identical to the one obtained by Voloshin 



4 Some Comments to the Moments 

In this section we will spend some time to discuss some interesting properties of the theoretical 
moments P^ h which have been calculated in Section |i| We will address three issues: (i) the relation 
between the strong dependence of the moments on Mj, and a s and the dependences of the moments 

On the SCaleS /i S oft; /-%ard 

and /Xfac, (ii) the properties of the resonance and continuum contributions 
and (iii) the quality of the nonrelativistic expansion. 

It is a characteristic feature of the moments that they depend very strongly on the bottom 
quark mass Mj, and the strong coupling a s . This is illustrated in Tab. || where the moments P^ 1 are 
displayed for n = 4,6,8,10,20 and for various values of Mj and a s (M z ) while the renormalization 
scales are fixed to /x so ft = 2.5 GeV and /ihard = Mfac = 5 GeV. The dependence on Mj is powerlike 
(P^ h ~ M fe _2n ) for dimensional reasons (see definition (0)). The dependence on a s is exponentially 
(see e.g. Eq. (|56|)) and comes from the resummations of the ladder diagrams containing the exchange 
of longitudinal Coulomb gluons. At this point one might conclude that fitting the theoretical moments 
to the experimental ones would allow for an extremely precise extraction of M& and a s , in particular if 
n is chosen very large. Unfortunately this conclusion is wrong. It is wrong from the conceptual point 
of view because for increasing n the effective smearing range AE in the integral (Q) becomes smaller 



27 



Moment 


Atsoft/[GeV] 


^hard/[GeV] 


Wac/[GeV] 




1.5 


2.5 


3.5 


2.5 


5.0 


10.0 


2.5 


5.0 


10.0 


Pf/[10- 8 GeV- 8 ] 


0.94 


0.37 


0.27 


0.31 


0.37 


0.43 


0.45 


0.37 


0.25 


Rf7[10- 12 GeV" 12 ] 


1.16 


0.41 


0.28 


0.34 


0.41 


0.47 


0.51 


0.41 


0.27 


Pf/[10- 16 GeV" 16 ] 


1.53 


0.49 


0.33 


0.41 


0.49 


0.56 


0.62 


0.49 


0.32 


P* V[10- 20 GeV" 20 ] 


2.10 


0.61 


0.39 


0.51 


0.61 


0.70 


0.79 


0.61 


0.39 


Ho /[ID" 40 GeV" 40 ] 


11.89 


2.37 


1.28 


1.98 


2.37 


2.72 


3.17 


2.37 


1.47 




Aihard = 5 GeV 
A*fac = 5 GeV 


/isoft = 2.5 GeV 
Mfac = 5 GeV 


//soft = 2.5 GeV 
/"hard = 5 GeV 



Table 3: The theoretical moments P^ h for n = 4,6,8,10,20 and fixed a s (M z ) = 0.118 and M b = 
4.8 GeV for various choices of the renormalization scales /u so ft, /"hard and /if ac . The two-loop running 
for the strong coupling has been employed. 



and smaller, which makes the perturbative calculations for the moments become less trustworthy [14]. 
In Section [2] we have used this argument to determine an upper bound on the allowed values on n. 
However, besides the conceptual arguments, the perturbative series for the moments itself contains 
a mechanism which prevents an arbitrarily precise determination of M b and a s for large values of n. 
In Tab. || the theoretical moments P^ h , n = 4,6, 8, 10, 20, are displayed for different choices for // so ft, 
/^hard an d /ifac an d for a s {M z ) = 0.118 and M b = 4.8 GeV. It is obvious that the dependence of the 
moments on the renormalization scales, and in particular on the soft scale, is becoming increasingly 
strong for larger values of n. As an example, the moment P^q (Pfo) can change by a factor of ten 
(five) if the soft scale is varied between 1.5 and 3.5 GeV. These huge scale dependences are mainly 
caused by the large NNLO contributions to the large n moments coming from the two-loop corrections 

(2) 

to the Coulomb potential, Vc , the second iteration of one-loop corrections to the Coulomb potential, 
V c {1) , and the non-Abelian potential, V^ A - During the fitting procedure, when all renormalization 
scales are scanned of the ranges (17), the large scale dependencies effectively compensate the strong 
dependence of the moments on M b and a s . In Section 6A it is shown that this affects mostly the 
extraction of a s rendering the sum rule, at least at the present stage, a rather powerless tool as 
far as precision determinations of the strong coupling are concerned. We want to stress that this 
compensation represents a very delicate balance which, if at all, can only be trusted if n is not chosen 
too large. We believe that this balance is still under controll for the values of n used in this work 
(4 < n < 10), although no proof for this assumption can be given. However, it is certain that for even 
larger values of n the extracted values of M b and a s might contain sizable systematic errors. 

We also would like to make one comment on the fact that the theoretical moments contain 
contributions from below (E < 0) and above (E > 0) the threshold point. As shown in Eq. (|32|), 
the former contributions come from the resonance poles whereas the latter arise from the continuum. 
To demonstrate the size of the resonance and the continuum contributions let us examine the LO 
contribution to Q n> \ which respect to this aspect. The contributions to [£> n ,i] LO from E < and E > 
can be calculated separately from Eq. (|48|) using the LO nonrelativistic expression for the cross section 
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0.0 


0.1 


0.2 


0.3 


0.4 


0.5 


0.6 


0.7 


0.8 


0.9 


1.0 


[f?n,l]B<0 


0.00 


0.02 


0.14 


0.50 


1.25 


2.65 


5.05 


9.01 


15.42 


25.66 


42.00 




1.00 


1.41 


1.92 


2.48 


3.09 


3.73 


4.39 


5.06 


5.74 


6.43 


7.13 



Table 4: The resonance (E < 0) and continuum (E > 0) contributions to the function [^ ni i] LO for 
0.0 < 6 < 1.0. 



from Eq. 
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In Tab. ||] expressions ( |68|) and ( |69[ ) are evaluated for 0.0 < </> < 1.0. For ~ 0.5 resonance 
and continuum contributions are approximately equal in size, whereas for larger values of <fi the 
resonance contributions dominate. This shows explicitly that for large values of n (where n > 4 can 
already considered as large) the resonance effects cannot be neglected. In particular, any sum rule 
analysis which is based on the large n moments and ignores the resonance contributions will lead 
to a bottom quark mass which is too low. From Eq. ( |68| ) it is also conspicuous that there are no 
resonance contributions proportional to a™ with n = 0, 1,2. This shows that in conventional multi- 
loop perturbation theory bound state contributions to the heavy-quark-antiquark production cross 
section (in lepton pair collisions) are produced by Feynman diagrams containing three and more loops. 
Because the twodoop level represents the current state of the art in covariant multidoop calculations 
where the full quark mass and energy dependence is taken into account (see Ref. [O] for a review 



and [42] far a recent publications on this subject), these peculiar contributions to the cross section 
sitting below the threshold point have not been observed so far. At the threedoop level, however, the 
cross section will have singular contributions oc a s /v 2 below the threshold. In fact, these contributions 
are required by the analyticity of the vacuum polarization function in the nonrelativistic regime. [Of 
course, for a proper description of the bound state regime fixed order multidoop perturbation theory 
is insufficient and a resummation of the singular terms to all orders in a s has to be carried out.] 

Finally, we also would like to address the question how well the nonrelativistic (and asymptotic) 
expansion at NNLO for the cross section R (Eq. (|5|)) and the integration measure ds/s n+l in the 
dispersion integral (f|) can approximate a complete covariant calculation of the large n moments, 
where all mass and energy dependences would be accounted for exactly. Strictly speaking, this question 
cannot be answered entirely because a complete covariant calculation of the moments, Eq. (fi), for large 
values of n is certainly an impossible task. (If it were possible, we would not use the nonrelativistic 
expansion and NRQCD in the first place.) However, a partial answer can be given by comparing the 
terms proportional to a™ with n = 0, 1, 2 in P„ , Eq. fl50D, to the corresponding contributions calculated 
in full QCD. For simplicity we only present a comparison of the Born and onedoop contributions in 
the following. The twodoop contributions lead to the same conclusions. The Born and the onedoop 
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contributions from Pi h read 
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The complete covariant versions of expressions ( |70|) and ( ffl| ) in full QCD can be determined from the 
well known Born and one-loop formulae for the cross section p0| , [H]], 
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where /? = (1 - AM^/q 2 ) 1 / 2 and p = (1 - /3)/(l + 0) and Li 2 is the dilogarithm, and the covariant 
form of the dispersion relation for the moments, Eq. (H), 



A 



71, QCD 



A 1 loop 

a fi,QCD 



ZN C Q\ 



7T 



4(4M 2 )™ra 3 / 2 



4M, 2 





3N C Q 2 ^ ,C F a 



A(4M^) n n 3 / 2 V vr 



s n+l 



1 -1 



P Born (s) 



s ra+l 



i? 



1 loop 



(«). 



(73) 



(74) 



4M2 



zSorn A 

,nrqcd; ^ 



A 1 loop 
n,QCD) T1,NRQCD 



Expressions (^) and ( |7"4| ) can be easily calculated numerically. In Tab. |B| A Bc 
and Ajj are presented for n = 1, . . . , 10. The difference for the Born (one-loop) contributions 
amounts to 6% (7%) for n = 4 and quickly decreases for larger values of n. Thus, for the values 
of n employed in this work the asymptotic expansion in the velocity and, in particular, the use of 
NRQCD, lead to a sufficiently good approximation to the exact covariant results for the cases where a 
comparison can be carried out. [At this point one has to compare the quality of the approximation to 
the large scale variations of the moments discussed at the beginning of this section.] This strengthens 
our confidence that our method to calculate the theoretical moments is sufficient at the level of 
the remaining theoretical uncertainties. In particular, we cannot confirm the claims in || that the 
nonrelativistic expansion would behave badly and would represent a good approximation only for 
n ~ 100. 
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n 


1 


2 


3 


4 


5 


6 


7 


8 


9 


10 


A Born 

n,NRQCD 


0.13 


0.56 


0.71 


0.78 


0.83 


0.85 


0.88 


0.89 


0.90 


0.91 


A Born 

a n,QCD 


0.60 


0.73 


0.79 


0.83 


0.86 


0.88 


0.89 


0.91 


0.91 


0.92 


A 1 loop 
^ll,NHQCD 


0.78 


4.25 


6.29 


7.87 


9.21 


10.41 


11.49 


12.50 


13.44 


14.33 


A 1 loop 

^n,QCD 


2.28 


4.25 


5.91 


7.36 


8.65 


9.84 


10.92 


11.94 


12.90 


13.81 



Table 5: The Born and one- loop contributions to the theoretical moments calculated in the nonrela- 
tivistic expansion (NRQCD) at NNLO and in full QCD for n = 1, . . . , 10. 



nS 


15 


25 


35 


M n5 /[GeV] 


9.460 


10.023 


10.355 


WlkeV] 


1.32 ±0.04 ±0.03 


0.52 ±0.03 ±0.01 


0.48 ±0.03 ±0.03 


nS 


45 


55 


65 


M n5 /[GeV] 


10.58 


10.87 


11.02 


r„ 5 /[keV] 


0.25 ±0.03 ±0.01 


0.31 ±0.05 ±0.07 


0.13±0.03± 0.03 


= q-^(10 GeV) = 131.8(1 ± 0.005) , (y/s) B § = 2 x 5.279 GeV 



Table 6: The Experimental number for the T masses and electronic decay widths used for the 
calculation of the experimental moments P^ ■ For the widths the first error is statistical and the 



second systematic. The errors for T15 and T25 are taken from [43]. All the other errors are estimated 



from the numbers presented in [44]. The errors in the T masses and the BB threshold (y/s) B § are 
neglected. 



5 Experimental Moments the Fitting Procedure 

In this section we will describe how the moments are calculated from experimental data and present 
our method to fit the experimental moments, P^ x ' , to the theoretical ones, P^ h . 

The experimental moments are determined using the available data on the T masses Mx( n 5) 
and electronic partial widths Fr(nS) = T(T(n5) — > e + e~) for n = 1, ... ,6. For a compilation of all 
experimental numbers see Tab. |6[ The formula for the experimental moments reads 

6 00 
9vr ^ T kS f ds 



+ / -^ttW(s). (75) 



n a 1 ^ M 2n+1 / s n+1 

V S BB 

The first term on the RHS of Eq. (fT5|) is obtained by using the narrow width approximation for all 
the known resonances 

Rres(s) = ^Y, T nS M nS 5( S - M* s ) . (76) 
Or 
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a em is the electromagnetic running coupling at the scale 10 GeV (see Tab. ||) which divides out the 
effects of the photonic vacuum polarization contained in the electromagnetic decay width. The 
second term describes the contribution from the continuum above the BB threshold. We approximate 
the continuum cross section by a constant with a 50% error 

rcont(s) = r c (1 ± 0.5) . (77) 

This simplifies the treatment of experimental errors in the continuum regime significantly but also 
represents an reasonably good approximation because for n < 4 the continuum is already sufficiently 
suppressed that a more detailed description of it is not necessary. During the fitting procedure we 
vary the constant r c between 0.5 and 1.5 which certainly covers all the experimental uncertainties. 
[In fact, this prescription renders the resonances 45, 5S and 6S, which lie above the BB threshold 
practically irrelevant.] 

For the fit we use the standard least squares method as described in The x 2 function 

which has to be minimized reads 

X 2 (M b ,a s )= £ {P t n h -P C n X ){S~ l ) nm {P t J:-P^)- (78) 
{n},{m} 

{n} represents the set of n's for which the fit shall be carried out and 5 _1 is the inverse covariance 
matrix describing the experimental errors and the correlation between the experimental moments. To 
construct the covariance matrix we use the errors in the electronic decay widths (where statistical and 
systematic errors are added quadratically) , in the electromagnetic coupling a em (see Tab. [fj) and the 
error in the continuum cross section, Eq. ([n|), which we also treat as experimental. The tiny errors 
in the T masses are neglected. At this point it is important to note that the errors in the electronic 
widths are certainly not uncorrelated due to common sources of systematic errors in the e + e~ collider 
experiments (mostly CLEO) where the widths have been determined. Unfortunately an analysis of 



these correlations cannot be found in the corresponding publications (see [44] for references). We 
therefore assume that the correlations between two widths can be written as 



5T nS 6T mS = acor ST s n y s s 5r s * s s , (79) 

where 5T n s is the systematic error in the electronic width T n s as given in Tab. ^ and a cor is a parameter 
which allow to switch the correlation on and off to check its impact on the extraction of Mb and a s . 
During the fitting procedure a cor is varied between zero (no correlation) and one (complete positive 
correlation of all systematic errors). Collecting all the quantities for which we take experimental errors 
into account into the vector 

in = { ris,r 2 s,r3s,r45,r 5 s,r 6 5,a em ,r COJ1 i j ; * = 1, — , 8 , (so) 

and using the standard error propagation formulae (see e.g. [£H|]) the covariance matrix reads 



° a pex 
q _ \ " ur n 

• • 1 OVA 



ftpex 



y °Vj 



Vij , (81) 



y 



10 To be more accurate, the electromagnetic coupling should be evaluated for each resonance individually at the 
corresponding resonance mass. The resulting differences, however, are smaller then the assumed error in a em itself and 
therefore neglected. 
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where 



/ (Wis) 2 ST 1S ST 2S 
5T 2S 5T 1S (5T 2S ) 2 



V; 



5F G s Wis 





6F G s ST 2 s 





W 1S 5T es 
W 1S 5T 6S 



(5a e 






{5r c f J 

\{j indicates that the functions are evaluated at the corresponding central values. 



(sr 2S ) 
o 










(82) 



The symbol 

The fitting procedure if complicated by the fact that the theoretical uncertainties (coming 
from the dependence of the theoretical moments on the renormalization scales « so ft, "hard and //f ac ) 
are much larger than the experimental errors, which are dominated by the errors in T\$, T 2 $ and T%s- 
Further, while it is reasonable to assume that the errors in the experimental data can be treated as 
Gaussian, this is certainly not the case for the "uncertainties" (or better "freedom") in the choices of 
the renormalization scales for which just a "reasonable" window can be given. It would therefore be 
inconsistent to include the theoretical uncertainties into the covariance matrix S. Nevertheless, it is 
important to have some means to combine both types of errors, the experimental and the theoretical 
ones. In this work this is realized by scanning all scales over the ranges given in Eqs. (|l7|). We will 
carry out two kind of fits. First, we fit for and a s simultaneously without taking into account any 
constraints on a s , i.e. ignoring all existing determinations of the strong coupling (Section |6.1|) , and, 
second, we fit for Mj, assuming that a s is a known parameter i.e. taking into account a constraint on 
a s (Section |6.2|) . 

To fit for Mb and a s simultaneously we employ a strategy closely related to the one suggested 
by Buras [45] and adopted by the BaBar collaboration |^6[ as a method to extract Cabibbo-Kobayashi- 
Maskawa matrix elements from various B-decays. Our strategy consists of the following two steps: 

(a) We first choose the range over which the renormalization scales u so ft, Uhard 

and //fac have to be 

scanned individually. For convenience we also count the constant r c , the correlation parameter 
a cor and the various sets of n's for which the fits shall be carried out as theoretical parameters. 
The individual ranges employed in this work are as follows, 



1.5 GeV 


< 


Msoft 


< 


3.5 GeV 


2.5 GeV 


< 


/^hard 


< 


10 GeV 


2.5 GeV 


< 


Mfac 


< 


10 GeV 


0.5 


< 


r c 


< 


1.5 





< 


(Icov 


< 


1. 



The sets of n's for which we perform the fits are 

{n} = {4, 5, 6, 7}, {7, 8, 9, 10}, {4, 6, 8, 10} 



(83) 



(84) 



The scanning over the ranges and sets given above is carried out by using a Monte-Carlo gener- 
ator. 
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(b) Then, for each set of theoretical parameters 



M = {// S oft , Mhard , Wac , r c , a cor , {n}} , (85) 

called a "model" , we construct the x 2 function as describe before and determine the 95% confi- 
dence level (CL) contour in the Mfe-o; s -plane by calculating the minimum x 2 -, X 2 m im an d drawing 
the contour x 2 (Mbi &s) = xLn + ^ The external envelope of the contours obtained for all models 
generated by scan represent the "overall 95% CL contour" which we will refer to as the "allowed 
range for and a s ". It should be mentioned that we do not impose a x 2 cu t which would 
eliminate models for which the probability of Xmin would be smaller than 5%. We will come 
back to this point in Section |(| 

We would like to emphasize that the allowed region for M& and ot s obtained by the procedure described 
above should not be understood in any statistical sense. In fact, it is quite difficult to ascribe any 
accurately defined meaning to the allowed region at all without reference to the method how it has 
been obtained. This is a consequence of the fact that the theoretical uncertainties, which cannot be 
apprehended statistically, dominate over the experimental ones. 

For the fit for Mb where a s is assumed to be a known parameter we treat a s like the theoretical 
parameters /x so ft , /^hard , A*fac , r c , a cor and {n}, i.e. we also scan over the given range of a s . The fit for 
Mb is then carried out in the same way as for the unconstraint fit described before. The only difference 
is that in this case the 95% confidence level "contour" for each model is determined by the equation 
X 2 (Mb) = Xmin + ^ because this method does represent only a one parameter fit. Some more remarks 



to this method can be found in Section 6.2 



6 Numerical Results and Discussion 



In this section we present the numerical results for the bottom quark pole mass M& gained from fitting 
the theoretical moments at NNLO calculated in Section || to the experimental moments obtained from 
experimental data. In Section 6T we discuss the result if M& and a s are fitted simultaneously ("un- 
constraint fit" ) and in Section 6^2 we present the result for Mj, if a s is taken as an input ( "constraint 
fit"). 



6.1 Determination of M& and a s without Constraints 

The result for the allowed region for and a s when both parameters are fitted simultaneously 
and no previous determination of a s is taken into account is displayed in Fig. [H]. The gray shaded 
region represents the allowed region in the M{,-a s plane. To illustrate that the allowed region does 
not have any well defined statistical meaning we have also shown the dots representing the best fits 
(i.e. the points in the Mb-a s plane with the lowest x 2 value for a large number of models. In fact, the 
region covered by the dots for the best fits is a measure for the size of the theoretical uncertainties 
inherent to our result. The latter uncertainties, which cannot be apprehended statistically, clearly 
dominate over the experimental ones, which are contained in the grey shaded region not covered by 
any dots. For convenience of the reader we have shown the result for a s at the scale \i = 2.5 GeV 
(lower frame axis) and \x = M z (upper frame axis) where we have used two-loop running for the strong 
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a s (M z ) 
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0.35 



Figure 10: Result for the allowed region in the Mb-a s plane for the unconstraint fit based on the 
theoretical moments at NNLO. The grey shaded region represents the allowed region. Experimental 
errors are included at the 95% CL level. The dots represent point of minimal x 2 f° r a large number 
of models. 



coupling. From the shape and orientation of the gray shaded region in Fig. [h] it is evident that Mb 
and a s are positively correlated. This can be easily understood from the fact that the theoretical 
moments are monotonically increasing functions of a s but monotonically decreasing functions of M& 
(see Tab. ||). However, we refrain from presenting a numerical value for the correlation because, as 
already mentioned, the allowed region for Mb and a s does not have any statistical meaning. 
For the bottom quark pole mass and the strong coupling we obtain 

4.74 GeV < M b < 4.87 GeV , (86) 

0.096 < a s {M z ) < 0.124, (87) 

0.175 < a s {2.5 GeV) < 0.308. (88) 
(NNLO analysis, Mb and a s are fitted simultaneously) 

Because the uncertainties for Mb and a s are not Gaussian we only present the allowed ranges obtained 
from Fig. |l(]. We would like to emphasize that in this context the inequality sign "<" does not have 
any mathematical meaning. It is only used to describe the bounds on M& and a s which are obtained 
from our fitting procedure. The allowed range for Mb, which spans over 120 MeV, can be definitely 
called a precise determination of the bottom quark pole mass. The allowed range obtained for the 
strong coupling, on the other hand, is consistent with the current world average, but much wider than 
the uncertainties of the latter. In addition, most of the allowed range for a s is located below the 
current world average. Taking the size of allowed ranges for Mb and a s as their uncertainty we arrive 
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at 

AM b 

M b 

Aa s (M z 



a s {M z ) 

Aa s (2.5 GeV) 
a, (2.5 GeV) 



2.5 % , (89) 
25 % , (90) 
50 % , (91) 



for the relative uncertainties in our determination of M b and a s . It is evident that the sum rule based 
on the large n moments, Eq. (||), is much more sensitive to the bottom quark mass than to the strong 
coupling. At least at the present stage one can certainly conclude that this sum rule does not belong 
to most powerful methods to determine a s as far as precision is concerned. 



From Eqs. ( |8q ) and (87) we can calculate the value for the running bottom quark mass. Using 



the two- loop relation between the pole and running mass |47j (see also Ref. [41| and references therein) 
and taking into account the correlation between the pole mass and the strong coupling we get 

4.09 GeV < m b (M r(ls) /2) < 4.32 GeV , (92) 
4.17 GeV < m b (m b ) < 4.35 GeV . (93) 

This result is in excellent agreement with a recent determination of the running bottom quark mass 
obtained from the three-jet rate in bb events at the CERN e + e~ LEP experiment DELPHI J48| , |48[ ], 
m&(Mx(ig)/2) = 4.16 ±0.14 GeV. The uncertainty in the result for the running quark mass, Eq. fl93|), 



is larger than for our pole mass result, Eq. (86), because of the correlation between M b and a s , which 
has to be taken into account in the conversion formula. 

We have checked that the allowed region for M b and a s presented in Fig. |l^ is insensitive to 
the particular choices of the scanning ranges for the renormalization scale //hard and the constants 
a cor and r c , which parameterize the correlation of the experimental data for the electronic widths and 
the continuum cross section above the BB threshold, respectively. However, the results depend on 
the choice of the ranges for the soft scale /i so f t and the factorization scale /if ac . This dependence is 



illustrated in Fig. 11 where we have displayed points for the best fits (a) for models with 1.5 GeV < 
Msoft < 2.5 GeV and 2.5 GeV < ^ soft < 3.5 GeV and (b) for models with 2.5 GeV < ^ fac < 5.0 GeV 
and 5.0 GeV < fi so ft < 10 GeV with different symbols. In both figures the other parameters have been 



scanned over the ranges given in Eqs. (J83|) . From Fig. 11(a) we see that the allowed range for M b 
does not depend significantly on the choice for the soft scale, whereas the allowed range for a s tends 
toward larger values if the soft scale is larger. Fig. |TT|(b), on the other hand, shows that the size of the 
allowed range for M b could be reduced if smaller factorization scales would be chosen. In that case 
the allowed range for a s would be only mildly affected. From this observation it might be tempting 
to choose the scanning range for /i so ft at higher scales and for /u.f ac at lower scales because this would 
lead to a seemingly more precise determination of M b and higher values for a s . However, we take 
the position that the choice of the scanning ranges for the renormalization scales should not depend 
on such considerations to represent a "reasonable choice". In fact, we consider it inappropriate to 
tune or "optimize" renormalization scales in some specific way if no good physical reason for that 
can be given. In our case the choice for the scanning ranges for the soft scale was motivated by the 
fact that it governs the nonrelativistic correlators for which (at NNLO) the relative momentum of the 
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Figure 11: Typical distribution of points representing the best fits (a) for models with 1.5 GeV < 
Msoft < 2.5 GeV and 2.5 GeV < /x soft < 3.5 GeV and (b) for models with 2.5 GeV < /i fac < 5.0 GeV 
and 5.0 GeV < /u so ft < 10 GeV based on the theoretical moments at NNLO. The other parameters are 
scanned over the ranges given in Eqs. (|83|). 



bottom quarks (which is of order M^a s ) represents the only relevant physical scale. Our choice for 
factorization scale /if ac , on the other hand, is inspired by the belief that is can take any value between 
the relative momentum of the bottom quarks and the hard scale which is of order the bottom quark 
mass (see Section |2|) . We will come back to this issue in Section |7[ 

It is very interesting to compare the results of our NNLO analysis presented above to an 
analogous analysis based on the NLO moments, i.e. ignoring all the NNLO contributions. [See the end 
of Section || for a prescription how the NLO moments can be recovered from the NNLO ones.] The 
result for the allowed range for Mb and a s based on the NLO moments is displayed in Fig. The 
gray shaded region and the dots have been obtained in the exactly the same way as described for the 
NNLO analysis. For comparison we have also indicated the allowed region obtained from the NNLO 
analysis by a polygon. Evidently the allowed region for Mb and a s covers a much larger area for the 
NLO analysis than for the NNLO one. At NLO the result for bottom quark pole mass and the strong 
coupling read 
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Figure 12: Result for the allowed region in the Mf,-a s plane for the unconstraint fit based on the 
theoretical moments at NLO. The grey shaded region represents the allowed region. Experimental 
errors are included at the 95% CL level. The dots represent point of minimal x 2 f° r a large number 
of models. The star and the diamond represent the results obtained by Voloshin and Kiihn et 
al. Q, respectively. The error-bars quoted by Voloshin are smaller than the symbol used to display 
his central value. The polygon represents the allowed region obtained from the NNLO analysis. 



4.64 GeV < M b < 4.92 GeV , (94) 

0.086 < a s (M z ) < 0.132, (95) 

0.144 < a s {2.5 GeV) < 0.368. (96) 
(NLO analysis, Mj, and a s are fitted simultaneously) 



From Fig, 12 and Eqs. ( |94] ) - ( J9q ) it is evident that the inclusion of the NNLO contributions of the 
moments leads to a considerable improvement upon a pure NLO analysis. We would like to point out 
that the uncertainties in Mj, and a s from our NLO analysis are much larger then the uncertainties 
quoted by Voloshin || and Kiihn et al. @. For comparison we have also displayed the results from 
Refs. ||, H] in Fig. 12. Because the theoretical moments used in Refs. (2|, ||] and the NLO moments 



used to generate the allowed region for Mj, and a s displayed in Fig. [T^ are equivalent, we consider the 
small uncertainties quoted in Refs. & Q as a consequence of an inappropriate treatments of the large 
theoretical uncertainties inherent to the perturbative calculations of the moments. (See Section ^ 
for a more detailed discussion.) Another way to see that that the NNLO contributions lead to a 
considerable improvement is to compare the distributions of best \ 2 values which are achieved by 
the models based on NNLO and NLO moments, respectively. In Tab. [7] the fraction (in percent) of 
best x 2 values within certain intervals is displayed for the NNLO and the NLO analysis based on, 
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NNLO 


28% 


17% 


16% 
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8% 


4% 


2% 


3% 


0% 


NLO 


0% 


0% 


0% 


0% 


0% 


1% 


7% 


35% 


57% 



Table 7: Distributions of best x 2 values for a NNLO and NLO analysis based on, at each case, 1300 
randomly generated models within the ranges (p3|). 



at each case, 1300 randomly generated models within the scanning ranges in Eqs. 83. Whereas for 
the NNLO analysis more than 60 % of the models have a best \ 2 value below 10, the bulk of the 
best x 2 values for the NLO analysis is larger than 50. We would like to emphasize that, because the 
uncertainties of the analysis are dominated by theory, the distributions of best x 2 values in Tab. 
represent only a measure for the quality of the theoretical expression for the moments, but do not 
contain any statistical information. We therefore cannot impose a x 2 on the models, let us say, based 
on an assumed statistical distribution of x 2 values. As an example, for two degrees of freedom and at 
the 95% CL, and assuming a Gaussian distribution such a x 2 cu t would eliminate all models whose 
best x 2 value is larger than 6. Evidently, in this case, none of the models based on the NLO moments 
would survive and we would be forced to reject, at least, the nonrelativistic expansion up to NLO as 
a legitimate tool to calculate the moments from QCD for the sets of n's considered in this work. 

6.2 Determination of Mb with Constraints on a s 

We now carry out the fitting procedure for Mb if a s is taken as an input, e.g. from the current world 
average. At this point one might be tempted to simply cut out of the gray shaded region in Fig. |l^ the 
part for which a s is located in the preferred range. Due to the sizable correlation between Mb and a s 
this would then lead to a much smaller uncertainty in Mb than given in Eq. (^6|). However, the naive 
procedure just described is not the correct way to account for a constraint on a s . This comes from the 
fact that for the unconstraint fit performed in Section |6.1| the strong coupling is essentially a function 
of the model parameters M. = {/U so f t , //hard> A'fao r c ^cor, {n}}, i.e. a s is not independent of the choice 
for M.. If a s is taken as an input, however, we have to treat a s and M. as independent, because we 
have to be able to freely assign values to them. Thus, if we take a s from the world average, we can 
expect that for a number of models the allowed range for Mb will be located outside the gray shaded 
region in Fig. [n]. As a consequence, the constraint fit will in general lead to larger uncertainties in 
Mb than the unconstraint one. In addition, due to the positive correlation between Mb and a s we can 
also expect that the result for the allowed region of Mb for the constraint fit will be located at slightly 
larger masses than for the unconstraint fit. 

We would like to point out that there are many ways to account for a constraint on a s which 
all might lead to slightly different results. In this work we account for a constraint on a s by treating it 
in the same way as the parameters in A4, i.e. we also scan over the preferred range of a s . The allowed 



range of Mb is then obtained in the same way as for the unconstraint fit carried out in Section 6.1 
with the difference that now only a one-parameter fit is performed (see also Section ||) . It should be 
noted that this method treats a s entirely on the same footing as the theoretical parameters in A4, 
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Figure 13: Result for the allowed M& values for a given value of a s . The grey shaded region 
corresponds to the allowed ranges for the NNLO analysis and the striped region for the NLO analysis. 
Experimental errors are included at the 95% CL level. It is illustrated how allowed range for Mj, at 
NNLO is obtained if 0.114 < a s (M z ) < 0.122 is taken as an input. 



i.e. the uncertainties on a s are not taken into account as Gaussian or statistical errors. In Fig. 13 the 
allowed range for Mb based on the NNLO moments is presented as a function of a s . For each given 
value for a s the allowed range for Mb, which is obtained by scanning all the parameters in A4 over the 
ranges (|83|), is the projection of the gray shaded region onto the Mb axis. If a region for a s is given 
the allowed range for Mb is obtained by projecting the gray shaded region for all the a s valued in the 
preferred region onto the Mb axis. As an example which is also illustrated in Fig. 13, staring from the 
world average for a s as given by Stirling Q, 0.114 < a s {M z ) < 0.122, we arrive at 



4.78 GeV < M b < 4.98 GeV (97) 
(NNLO analysis, a s (M z ) taken from the world average Q) 

for the bottom quark pole mass. The result is consistent with Eq. (^) obtained from the unconstraint 
fit. However, as expected, the allowed range for Mb is wider and, in addition, located at slightly larger 
masses. In fact, the uncertainty on M& for the constraint for is almost a factor of two larger. We 
have checked that the result for Mb is insensitive to the particular choice of the scanning ranges for 
^hard> /Afac> a cor and r c . However, the bottom quark mass tends toward lower values if // so f t is chosen 
larger. We have also displayed the result for the NLO analysis in Fig. ^ as the striped area. As 
for the unconstraint fit the inclusion of the NNLO contributions to the moments leads to a smaller 
uncertainty for Mb, although the improvement is not as dramatic. We want to mention that the larger 
uncertainty for Mb obtained from the constraint fit is partly a consequence of the fact that our fitting 
procedure does not treat the error on a s as Gaussian or statistical. Therefore one might argue that the 
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uncertainties in Eq. (^) are too conservative. However, from the way how a world average is gained, 
it is certain that the error of a s constrains a sizable systematic contribution. Because an accurate 
quantitative description of such a systematic error is quite difficult, we take the position that the error 
on a s should be treated in a conservative way. 

Using the result in Eq. (97) and the two-loop relation between the running and the pole 
mass [47] we obtain 

4.08 GeV < m b (M r{ls) /2) < 4.28 GeV, (98) 
4.16 GeV < m b (m b ) < 4.33 GeV. (99) 

for the running bottom quark mass. It is remarkable that this result and the result for the running 
quark mass based on the unconstraint fit, Eq. (|93[), are almost identical. 



7 Comments on Previous Analyses 

In the past few years there have been three previous analyses by Voloshin ||, Jamin and Pich || and 
Kiihn, Penin and Pivovarov Q where the bottom quark pole mass and the strong coupling have been 
extracted from data on the T mesons and using the same sum rule as in our analysis. We would like 
to emphasize that in Refs. ]2j, |3], ||] no consistent determination of NNLO corrections has been carried 
out and that the results by Voloshin 

M b = 4.827 ± 0.007 GeV , a s (M z ) = 0.109 ± 0.001 , (Voloshin) (100) 

Jamin and Pich (JP) 

M b = 4.60 ± 0.02 GeV , a s (M z ) = 0.119 ± 0.008 , (Jamin, Pich) (101) 

and Kiihn, Penin and Pivovarov (KPP) 

M b = 4.75 ± 0.04 GeV a s (M z ) = 0.118 ± 0.006 , (Kiihn et al.) (102) 

are contradictory to each other and partly to our own results. In particular, although no NNLO 
contributions have been included, all results in Refs. [|, [|] are claimed to have much smaller uncer- 
tainties than any of the results obtained in our analyses. In this section we will explain the origin of 
those discrepancies and give some comments on the methods used in Refs. [||, ||, |j] from the point of 
view of the strategies followed in this work. To organize the discussion we will analyze the methods 
used in Refs. [§, ||, ||| with respect to three aspects: (i) theoretical expression for the moments, (ii) 
optimization and tuning of the perturbative series for the moments and (hi) fitting procedure and 
error analysis. Because the theoretical uncertainties in the determination of M b and a s are much 
larger than the experimental ones, we will neither focus on the treatment of experimental errors nor 
on the formulae used for the experimental moments. Compared to the effects caused by using different 
methods to handle the theoretical uncertainties, the differences in the treatment of the experimental 
side of the analysis represent only a minor issue. We also would like to mention that in the analyses 
of Voloshin, JP and KPP moments with n as large as 20 were used. According to the estimates given 
in Section |2| this means that the effective smearing range contained in those moments is already of the 
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same size as or even smaller than A QCD . This leads to a an additional source of systematic theoretical 
errors in the results of Voloshin, JP and KPP. We have checked, however, that using moments with 
10 < n < 20 causes only shifts in the results for M& (and a s ) which are small compared to the size of 
the theoretical uncertainties at the NLO level as we have estimated them in our analysis. Therefore 
we will not raise this issue in the following discussion. For a NNLO analysis, however, where the 
uncertainties in M& are shown to be smaller, the use of values of n which are too large is an important 
issue and can lead to considerable errors. 

Theoretical expressions for the moments: 

Voloshin's moments are identical to ours at the NLO level. 

The moments used by KPP have been calculated in the same way as Voloshin's (and ours at the NLO 
level) with the difference that the dispersion relation in Eq. (||) has been performed numerically in 
terms of it covariant form, i.e. without using the asymptotic expansion (p7| ) and the inverse Laplace 
transform. We have checked that for the values of n's considered by Voloshin, KPP and us the dif- 
ference between both approaches is negligible. Thus, the moments used by KPP are equivalent to 
Voloshin's and ours at the NLO level. 

The moments by JP, on the other hand, were obtained from the Born, one-loop and two-loop ex- 
pressions for R{e + e~ — * bb) supplemented by a resummation of LO Coulomb singularities in form 
of the Sommerfeld factor (see Eq. (|3^)). Further, the one-loop corrections to the Coulomb potential 
have been implemented by inserting them directly into the Sommerfeld factor, i.e. without using time- 
independent perturbation theory. For the dispersion integration @ JP have only taken into account 
cm. energies above the threshold point (s > 4M?). We disagree with the moments used by JP in 
two major points. Most important, JP did not take into account the bound state poles of the cross 
section R, which are located below the threshold point (s < 4M^ ). We have demonstrated in Section || 
that the bound states represent the dominant contribution to the moments for large values of n (see 
Tab. Thus the moments used by JP are far too small which causes the bottom quark pole mass 
obtained from the fits to be too low.[^] In fact, one can easily see that omitting the bound state poles 
for large values of n will always lead to a bottom quark pole mass M& < Mwis)/2 4.7 GeV regard- 
less whether a s is determined from the fit or taken as an input. This explains why the value for M& in 
the analysis by JP is significantly smaller than in the analyses by Voloshin, KPP and us. In addition, 
we do not think that the effects of the running of the Coupling governing the Coulomb potential have 
been treated properly. JP simply inserted the one-loop corrections to the Coulomb potential into 
the Sommerfeld factor. Whereas this legitimate for the non-logarithmic corrections, it is not for the 
logarithmic ones because the effects arising from virtual momenta below and above the scale ~ Mba s 
are not taken into account correctly. This can only be achieved by using time-independent perturba- 
tion theory (or by solving the Schrodinger equation exactly). We therefore conclude that the results 
obtained by JP contain large systematic theoretical errors which are by far larger than indicated by 
their error analysis. That the value for a s obtained by JP still seems reasonable is a consequence of 
the fact that the moments are much less sensitive to a s than to 

Optimization and Tuning: 

n The same conclusion has been drawn in Ref. Ml . 
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We have shown in Section || that the perturbative corrections to the moments and the resulting scale 
dependences are quite large. This behavior is particularly obvious at the NNLO level. However, al- 
ready at NLO the corrections are uncomfortably large. In our analysis this feature has been fully taken 
into account during our fitting procedure. In fact, it is the main source of theoretical uncertainties 
in our results. In the analyses by Voloshin and KPP, however, the perturbative expansion for the 
theoretical moments has been tuned to improve the convergence. 

In Voloshin's work, at each value of n the soft scale //soft has been fixed such that the NLO corrections 
caused by Vc, Eq. ([TS|), vanish exactly and the hard scale /ihard has been fixed to the BLM scale 
Thus, Voloshin has eliminated the scale dependences of the moments. We would like to emphasize 
that we consider Voloshin's prescription as one possible choice for the renormalization scales, which 
essentially corresponds to selecting one single model out of the range of models used in our analysis. 



We have shown in Section y (see e.g. Fig. |ToD that the results for M& and a s depend significantly 
on such a choice. Because we think that no argument can be found why Voloshin's choice should be 
better than others, we have the position that a scan over all "reasonable" models should be carried 
out. Because Voloshin has not carried out such a scan we consider the theoretical uncertainties quoted 
in his analysis as largely underestimated. 

In the analysis by KPP, at each value of n a non-logarithmic piece of Vc has been absorbed into the 



LO nonrelativistic Green function, Eq. (31), such that the NLO corrections caused by the non-absorbed 
piece (calculated via first order time-independent perturbation theory) vanish. This optimization is 
quite similar to Voloshin's but leaves the soft scale unfixed. It should be mentioned that KPP have 
explicitly identified soft and hard scale which has eliminated the possibility to vary both scales inde- 
pendently. This reveals why the uncertainties quoted by KPP are much larger than Voloshin's, and 
partly explains why they are still much smaller than the uncertainties obtained from our NLO analysis 



where no optimization has been performed. (See Fig. 12 for a graphical comparison). 

JP have not carried out any optimization. However, due to their way to calculate the moments start- 
ing from the expressions of the covariant multi-loop expressions for the cross section, JP implicitly 
identified soft and hard scale. 

Fitting procedure and error estimate: 

In the analysis by Voloshin and KPP a two parameter fit was carried out to obtain Mf, and a s for 
the sets {re} = {8,12,16,20} and {10,12,14,16,18,20}, respectively. Thus, the results obtained by 



Voloshin and KPP should be compared with the results of our unconstraint fit presented in Section 6.1. 
Because Voloshin has eliminated all scale dependences, he has estimated the theoretical uncertainties 
in his analysis using the assumption that the NNLO corrections to the re'the moments can be parame- 
terized by a global factor (1 + c/n) where c is a number of order one. The size of the uncertainties was 
obtained from the variation of the best fits for Mb and a s if c is first fixed to zero and then obtained 
from a three parameter fit. The theoretical uncertainties gained by this method have been of the same 
size as the (small) experimental errors. We have shown in Section |I] that the NNLO contributions to 
the moments have an entirely different structure (large size, growing with re, tremendous dependence 
on the soft scale) and cannot be accounted for by the global factor (1 + c/re). Thus, Voloshin's method 
to estimate the theoretical error is not capable to account for the true size of theoretical uncertainties 
inherent to the perturbative calculation of the moments. 
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In the analysis by JP the theoretical uncertainties for Mb and a s were essentially obtained from the 
variation of the best fit for Mb and a s (for fixed /4 ft = Hh&rd = Mb) when the two-loop corrections 

(2) 

to the Coulomb potential, Vc , are included and when the two-loop contributions to the high energy 
cross section are removed. No additional uncertainties (e.g. from the renormalization scale depen- 
dence) have been taken into account based on the argument that this would lead to a double counting 
of theoretical uncertainties. We disagree with this statement because the effects of the inclusion or 
removal of the two-loop corrections to the high energy cross section or the Coulomb potential certainly 
depends on the value of the other parameters (like the renormalization scales). This and the fact that 
JP have neglected the bound state poles, which are the dominant source of large corrections to the 
moments (and their scale dependence) for large values of n, have lead to an underestimate of the 
theoretical uncertainties (besides the large systematic errors mentioned above). 

KPP, finally, have determined Mb and a s separately. For the determination of M& a s (M z ) = 0.118 
was taken as a fixed input. Thus, the result for M& by KPP should be compared to the results of our 
constraint fit presented in Section 3.2. The method used by KPP to obtain Mb was based on solving 
the equation P*" = P% x for Mb while n and all the other parameters are fixed to specific values. The 
mean value and the uncertainty for Mb has then been gained by calculating the mean and observ- 
ing the spread of Mb values when this procedure was carried out, first, for fi so f t = //hard = Mb and 
n = 10, 12, . . . , 20 and, second, for fixed n = 14 and 1.2 GeV < /i so f t = //hard < Mb- This procedure 
effectively scans over some fraction of the range of models used in our fitting procedure but misses 
e.g. models with n ^ 14 and /U so ft < Mb- This the main reason why the uncertainties quoted by KPP 
are much smaller than in our NLO analysis. 



From the discussion presented above we come to the following final conclusion about the results 
by Voloshin, JP and KPP in comparison to our own analysis: The theoretical moments calculated by 
Voloshin and KPP are equivalent to our NLO moments. We therefore consider the results obtained by 
Voloshin and KPP consistent with our own results at the NLO level (see Fig. ^ and 13). However, the 
theoretical uncertainties are underestimated in both analyses, which leads to the apparent contradic- 
tion between the results by Voloshin and KPP. In view of the error analysis performed in our analysis, 
where we tried to impose as less bias as possible, the results by Voloshin and KPP are perfectly 
consistent. The apparent contradiction essentially corresponds to a disjunct (and from our point of 
view biased) choice of models used for the fitting procedure in both analyses. We want to emphasize 
that the choice by Voloshin is not less plausible than the one made by KPP which illustrates the need 
that the whole range of models must be scanned in the fitting procedure. The comparison between 
the results of our analysis, where such a complete scan has been performed, and the results obtained 
by Voloshin and KPP makes this obvious. The theoretical moments determined by JP, on the other 
hand, do not take into account the bound state poles which represent the dominant contributions to 
the moments for large values of n. As a consequence the Mb value obtained by JP is too small and 
has to be considered inconsistent with the results by Voloshin, KPP and us and, in particular, with 
the nonrelativistic expansion of QCD, where the bound state poles are predicted. That the value for 
a s obtained by JP still seems reasonable is a consequence of the fact that the moments are much less 
sensitive to a s than to Mb- 

While this paper was in its final stages there appeared a letter by Penin and Pivovarov (PP) fo"l[j 
where the NNLO corrections to the large n moments have also been included in the sum rule deter- 
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mination of the bottom quark pole mass. The formulae for the moments used by PP were based on 
previous results for the cross section published in [Q, || [| 1C ] and are therefore conceptually equivalent 



to ours. For the bottom quark pole mass PP quote the result = 4.78 ± 0.04 GeV. The result is 
consistent with ours. The uncertainty, however, is smaller and the allowed range for M& is somewhat 
lower than for our NNLO results. To obtain this apparently more precise result PP have used the 
same methods as in Ref. [Bj (which we have already criticized above) with the difference that all the 
scales (including the factorization scale) were varied in the range ± 1 GeV. We consider this range 
too high for the soft scale. Further, there is a sign error in the CaCf piece of the O(a^) short-distance 
coefficient in Eq. (2) of [51], which is not contained in the original publication Q, and in Eqs. (42) and 
(ff5|) of this paper. These issues and the fact the PP used values of n for the moments between 10 and 
20, which we consider too large for a NNLO analysis, are the main reasons why the result by PP is 
located at lower masses. Further, we would like to point out two mistakes in |[l| which, however, do 
not affect the results for M^. PP state that the knowledge of the nonrelativistic correlators is sufficient 
to completely determine the vacuum polarization function, Eq. (|l|), in the threshold regime. This is 
not true due to UV-divergences in the real parts of the correlators which still have to be removed by 
an additional matching calculation. (See the discussion at the end of Section |3.2| .) PP state that the 
factorization scale //f ac in the 0(a 2 s ) short-distance constant and in the correlators (Eqs. (2), (5) and 
(6) in |)l|])) is defined in the MS scheme. This statement is not true. The formulae used by PP were 
taken from Ref. || where the factorization scale is defined in a cutoff scheme. The corresponding 
results in the MS scheme can be obtained by an additional redefinition of the factorization scale //f ac - 
(See the discussion following Eq. d3l"|).) 



8 Conclusions and Outlook 

Based on the argument of global duality and causality one can relate the derivatives of the vacuum 
correlator of two bottom-antibottom vector currents at zero momentum transfer to an integral over 
the total cross section of the production of hadrons containing a bottom and an antibottom quark in 
electron-positron annihilation 



n\ q 2 ( ill/' 2 ^ " 



% > ( 103 ) 



It is therefore possible to relate a theoretical calculation of the moments P n to experimental data for 
the cross section R(e + e~ — > "66 hadrons"). The limit of large values of n is of special interest for 
relation ( |103| ) because in this limit the high energy contributions are suppressed. For the theoretical 
side this means that the bottom-antibottom pair can be treated nonrelativistically, and for the exper- 
imental side that data for the production of T mesons are already sufficient to saturate relation ( |103| ). 



The requirement that the effective range of integration is larger than A QCD 14] imposes an upper 
bound on the values of n for which a perturbative calculation of the moments can be trusted. Due 
to the large size of the bottom quark mass of order five GeV we are in the fortunate situation that a 
range of values of n can be found for which the 66 system can be treated nonrelativistically and, at 
the same time, the range of integration is still broad enough compared to A QCD . We have identified 
this range of "large values of n" as 4 < n < 10. In this work we have used the arguments just given 
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to determine the bottom quark pole mass M& and the strong coupling a s in the MS scheme from 
experimental data on the T masses and electronic decay widths. 

The aim of this paper was twofold: 
1.) Calculation of NNLO corrections: The complete set of NNLO corrections in relation ( [103[) for 
large values of n has been calculated. This includes corrections to the expressions in the nonrelativistic 
limit or order a 2 , a s v and v 2 , where v is the velocity of the bottom quarks in the cm. frame. The 
conceptual difficulty in those calculations is that the relativistic corrections, e.g. from the kinetic energy 
or from higher order interactions like the Darwin or the Breit-Fermi potential lead to ultra-violet 
divergent integrations. We have used the concept of effective field theories formulated in NRQCD [11, 



12] to deal with this problem. In NRQCD the latter divergences appear as a natural consequence of 



the existence of higher dimensional operators which lead to the renormalization of lower dimensional 
ones. The exact form of the renormalization constants is obtained through matching to full QCD. 
This automatically provides a separation of all relevant effects into short-distance (contained in the 
renormalization constants) and long-distance ones (contained in matrix elements). In our case this 
leads to an expression for the moments (and the cross section R(e + e~ — ► "66 hadrons") which is a sum 
of terms each of which consists of a nonrelativistic current correlator multiplied with a short-distance 
factor (see Eqs. ( |l3|) and (|50|)). For the leading term in this series we have performed matching at 
the two-loop level. Although the NNLO contributions are quite large, they lead to a considerable 
reduction of the theoretical uncertainties in the extraction of Mb. 

2.) Conservative approach for the error estimates: The uncertainties in the determination of 
Mb (and a s ) based on sum rule ( |103| ) are dominated by theory, in particular, by the remaining large 
renormalization scale dependences of the theoretical moments. These scale dependences are caused by 
large coefficients which arise in the perturbative calculation of the corrections to the moments in the 
nonrelativistic regime. In contrast to statistical errors, which can be treated in a standardized way, it is 
not an easy task to properly account for theoretical uncertainties, in particular, if a model-independent 
(in the framework of QCD) analysis is intended. In fact, in an analysis where theoretical uncertainties 
dominate results may easily become biased depending on personal preferences. For the case of the 
determination of Mj, from sum rule ( |103|) this has lead to the paradoxical situation that in two earlier 
publications [^, |j] contradictory results were obtained although equivalent theoretical expressions for 
the moments were used. In this paper it was attempted to include as less personal preference into the 
analysis as possible by scanning all theoretical parameters independently over reasonably large ranges 
which were in size and location motivated from general considerations. For each set of parameters, 
called a "model", a standard statistical fitting procedure was carried out using the method of least 
squares to calculate a 95% CL contour. The external envelope of the contours obtained for all the 
scanned models was then taken as the "overall allowed range" which, we want to emphasize, does not 
have any well defined statistical meaning due to the dominance of theoretical uncertainties. This makes 
the scanning method more conservative (and in our opinion also more honest) than the methods used 
in Refs. [§, |J. In addition, the scanning method has the advantage that it automatically accounts for 
non-linear dependences on the theoretical parameters and prevents by construction the Gaussian-like 
treatment of theoretical uncertainties. Of course, the results presented in this work are not completely 
free from personal preferences either because of the choice of the ranges used for the scanning. 

In this paper we have performed two different analyses based on the scanning method and 
including the new NNLO corrections. First, and a s were determined simultaneously using the 
least squares method for two parameters. We have obtained 
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4.74 GeV < M b < 4.87 GeV , 

4.09 GeV < m b (M r(ls) /2) < 4.32 GeV, 
0.096 < a s (M z ) < 0.124 

(NNLO analysis, M b and a s are fitted simultaneously). 



(104) 
(105) 
(106) 



The corresponding result using the NLO expressions for the moments yielded considerably larger un 

and (P 



certainties (see Figs. |l0| and |l2| and Eqs. (|8q)-(|88|) and (|94[) ~(Pq)). The results show that relation ( |103[ ) 
allows for a much more precise determination of the bottom quark mass than for the strong coupling. 
Second, M b was determined using the least squares method for one parameter and taking a s as a 
known parameter. We have obtained 



4.78 GeV < M b < 4.98 GeV , (107) 

4.08 GeV < m b (M T(ls) /2) < 4.28 GeV (108) 

(NNLO analysis, 0.114 < a s {M z ) < 0.122). 

As for the first analysis the NNLO contributions to the theoretical moments lead to a reduction of 
the uncertainties (see Fig. |l^). In our opinion, the sum rule ( |103j ) can be regarded as a quite precise 
tool to determine the bottom quark mass. For the determination of the strong coupling we are by far 
less optimistic. 

In the past few years there have been three previous analyses by Voloshin 0] , Jamin and Pich ||] 
and Kiihn et al. Q where the bottom quark pole mass has been determined from experimental data 
for the masses and the electronic decay widths for the T mesons using the sum rule (|103|) . The results 
obtained in those three analyses are contradictory to each other, and, although NNLO corrections 
have essentially not been included, quote uncertainties smaller than in our own NNLO analysis. The 
results obtained by Voloshin and Kiihn et al. are based on moments which are equivalent to ours 
at the NLO level. In view of the uncertainties for M b (and a s ) obtained from our analysis at NLO, 
the results by Voloshin and Kiihn et al. can therefore be regarded as consistent with each other (and 
us), see Fig. [l^. The small uncertainties quoted by Voloshin and Kiihn et al. come from too tight, 
model-like bounds imposed on the theoretical parameters. The results obtained by Jamin and Pich, on 
the other hand, contain a large systematic error due to the negligence of the bound state contributions 
in the moments. We consider the result by Jamin and Pich inconsistent with those by Voloshin, Kiihn 
et al. and us, and in particular with the nonrelativistic expansion of QCD. 

It is quite interesting to ask whether and how the results determined in this work can be fur- 
ther improved in order to arrive at even smaller uncertainties for the bottom quark pole mass or the 
strong coupling. From the technical point of view the answer would simply by to calculate the NNNLO 



contributions in relation (103). Such a task, however, is highly nontrivial. Apart from the fact that a 
three-loop matching would have to be performed also the NNNLO effects in the bottom-antibottom 
interactions would have to considered. This would require a consistent treatment of retardation ef- 
fects which are caused by the non-instantaneous exchange of gluons and, as a prerequisite, a better 
understanding of higher order Fock bottom-antibottom-gluon states. In principle a calculation to 
determine these effects would be the QCD analogue of the determination of the Lamb shift contri- 
butions to the positronium wave function. So far no technical instruments have been developed yet 
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to immediately tackle this challenging problem. We further believe that it is unlikely that this goal 
can be achieved entirely in the framework of perturbation theory because it involves also the bound 
state energy ~ M& v 2 ~ Mj, a 2 , as a relevant scale. For the bottom quark this scale is already of the 
same size as the typical hadronization scale Aq CD , which means that the bottom-antibottom-gluon 
propagation is certainly nonperturbative. In fact, the rather uncomfortably large NNLO corrections 
in relation ( 103| ) might be regarded as a first warning sign in support of this view. 

At this point it seems to be just natural to mention the renormalon ambiguities contained in 
the definition of the pole mass [52] which is defined perturbatively as the location of the singularity 
of the renormalized quark propagator. This ambiguity indicates that the pole mass has an intrinsic 
uncertainty of order A QCD ~ 200 — 300 MeV. It is caused by the long range sensitivity of the pole mass 
and reflected in a factorial growth of the high order coefficients of the perturbation series connecting 
the pole mass to other mass definitions like MS which seem to be free from this problem. Our results 
and the rather pessimistic prospect to further improve the results obtained in this paper certainly 
support this view. However, the notion of the renormalons might also give hints toward a more precise 
determination of the bottom quark mass because it implies that with a different mass definition the 
perturbative series for the moments might become better behaved. In this work we have not attempted 
to make use of this possibility, but we hope to return to this issue in the near future. 
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A NNLO Corrections from 5Hkm, V BF and V NA 

In this appendix we present some details about the calculation of the NNLO corrections to the zero- 
distance Green function coming from the kinetic energy 5Hki n (r) = — V 4 /4M^, the Breit-Fermi po- 
tential Vbfj Eq. (|l8|), and the non-Abelian potential V NA , Eq. (|2l|). 

At NNLO the corrections coming from 5Hki n , V BF and V NA are determined from first order 
time-independent perturbation theory, 

[G f M «« » = - jSr^r,E )W ?Hr, a ,E), (109) 

where 

5H(r) = -V + V BF (f) + V NA (r) . (110) 

b 

Because the zero-distance Green function only describes bottom-antibottom pairs in a 3 S"i triplet state, 
we can take the angular average and evaluate the spin operators for 5H in expression ( |109| ). The form 
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of 5H then simplifies to 



V 1 



4M b 3 



C F a s V 2 11 C F a s ir (3) 
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2M b r 2 



(111) 



Using the equation of motion for the Coulomb Green function, Eq. fl2q), we can eliminate the V 
terms in SH^si- For illustration, let us consider the corrections coming from the term —Cp^j^z in 

b 

SHzsi ■ Using the equation of motion we arrive at the relation 
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The third term in the brackets represents a power divergence which is dropped in our convention (see 
the text after Eq. (j30|)). Using the same arguments for the kinetic energy term we get 
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Collecting all terms from Eqs. ( |lll|) -( |lT3|) we arrive at 

[G( 2 )(o,o,£)] kin+BF+NA 
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The first and the second term in the brackets on the RHS of Eq. ( |114j) are handled by redefining the 
energy, E —* + S 2 /4M 2 and the coupling, a s — > a s [l + 3£'/2Mft], in the nonrelativistic Coulomb 
Green function. The calculation of the 5-function term is trivial. The treatment of the 1/r 2 term, 
on the other hand, is rather awkward. However, we can infer the correction caused the 1/r 2 term by 
using the facts that the wave functions to the Schrodinger equation 



V 2 

m 2 



a 
r 



E *(r) = 



(115) 



can be determined exactly for any energy E (see e.g. Ref. j53|) and that the imaginary part of the 
Green function G of Eq. (115) in the continuum, i.e. for any positive energy E, is proportional to the 
modulus square of the scattering wave function at the energy E. From this it is straightforward to 
derive for positive energies the relation 



ImG(0,0,£) 



lim 



(2pr) 



mp 

4-7T 



cxp 



airm 
2p 



r(i + a -i|s 



T(2s + 2) 



(116) 
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where s(s + 1) = —b and p = y/m(E + Tej. Expanding the RHS of Eq. (|116| ) in small bf^\ and imposing 
the short-distance cutoff /if ac as described in Section |34] we obtain for positive energies the relation 



Im 



d 3 rG^(0,r, E) 



r 2 n 2 

( - / F a s 
M b r 2 



Gf\r,0,E) 



4 Cf cl s tt 



Im 



(g^(o,o,e) 



(117) 



Due to analyticity relation (|117| ) is then also valid for any real energy. Up to (irrelevant) constants 
we can therefore write 



d 3 fG^(0,r, E) 



r 2 n 2 
M b r 2 



Gf\r,0,E) 



4 Cf a s tt 
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(118) 



Collecting all terms the final result for the sum of the zero-distance Coulomb Green function and the 
NNLO corrections caused by SH^, V BF and V^ A reads 



G^ (0,0, E,a s ) + (0,0, E)]"* 



+BF+NA 



1 + 



E 
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+-. 



1 + 



3 Ca \ Cf o s " 
1 2C F ~ 
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G^ (0,0, E,a t 



2M b 

2 



(119) 



up to corrections beyond the NNLO level. g£° (0, 0, E, a s ) is defined as the expression on the RHS 



of Eq. (|3l|). Rewriting the energy in terms of v = y/(E + ie)/Mb we arrive at the result displayed in 
Eq. ©. 

We do not want to leave unmentioned that for the treatment of the singular 1/r 2 potential we 
have ignored the fact that its coefficient (mainly through the large non-Abelian contribution) is large 
enough that the bb system can collapse to a point (see e.g. p3| ). This would lead to the breakdown 
of hermiticity. Thus, the result in Eq. ( |119| ) has some heuristic character. However, we strictly treat 
the singular 1/r 2 (and also the 5^(r)) potential as a "small" perturbation to the Coulomb exchange 
and remove the arising UV singularities through the matching procedure. No exact treatment of the 
singular potential is intended. In this sense the result in Eq. (|119[) should be fine. 



B Inverse Laplace Transforms 

In this appendix we present the list of inverse Laplace transforms used to the calculate the theoretical 
moments at NNLO. In the following we use the conventions 

dlnT(z) 

V(z) = 



dz 

d n 

*(»)(*) = -=-*(*), 

w dz n v ' 



12 Because we want to treat the 1/r 2 potential as a perturbation, the limit r — > has to be taken after the expansion 
in b. 
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)F 2 (a,b;z) = r(o)T(6) £ 



- o r(a + A;)r(6 + A:) fc! 



7+ioo 
1 f 1 



27ri 7 a;" 

7— ioo 
7+ioo 

1 f Inx 



e^dx 



27ri 7 x 1 ' 

7— ioo 



e x dx 



7+ioo 

1 f In 1 



27TZ J .)•' 

7— ioo 



7+ioo 



2vri 



ln d x 



7—400 



7+ioo 



- 1 - / — sm(-?=)e xt dx 
2iri J x v K yfx' 



7— JOO 

7+ioo 



j r \w 3} d 

/ sm(—=)e xt dx 

2iri J x u y/x 



7—400 



7+ioo 



1 [ In x a t 
2iri 7 x v \ x 



7—100 



tV-l 



Y{y) 



e xt dx = ^— l - { -lnt] 2 - } 



$(1/) - hit 



r(» 



i> + 



— oF2( 2'^ 1 2' 



1 a 2 



at""* 



(-5) 



hxt 



, , 3 1 a 2 
oF2(^f+2»-jt 



d „ /3 la 2 



(120) 



(121) 



(122) 



e xt dx = f^{[^(^)- lnt ] - 3 In*] ¥» + ¥"(»/)}, (123) 
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C The Constants w^' 1,2 and w^ 1,2 

In this appendix the constants iv®' 1,2 and tip' 1,2 from expression (p7|) are given. They generically 
parameterize the higher order contributions to the Green function of the Schrodinger equation (24) 
coming from the radiative corrections to the Coulomb potential, Vc and vj 2 \ Eqs. (|19|) and (22). For 
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the constants Wp' 1,2 we were able to calculate analytic expressions. The results read (p = 1, 2, 3, . . .) 
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The constants ufy 1,2 are calculated numerically. The corresponding integrals are (i = 0, 1,2) 
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where 
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and 



x = 1 + t + u , 
y = 1 + w — s . 
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